Skip to content

Commit b0acbfd

Browse files
cgcgcgcgcgcg
committed
updates
1 parent 06c46ce commit b0acbfd

File tree

27 files changed

+1287
-399
lines changed

27 files changed

+1287
-399
lines changed

‎base/PyNucleus_base/tupleDict.pyx

+2-1
Original file line numberDiff line numberDiff line change
@@ -208,7 +208,8 @@ cdef class arrayIndexSet(indexSet):
208208
for i in s:
209209
self.indexArray[k] = i
210210
k += 1
211-
qsort(&self.indexArray[0], self.indexArray.shape[0], sizeof(INDEX_t), compareIndices)
211+
if k > 1:
212+
qsort(&self.indexArray[0], self.indexArray.shape[0], sizeof(INDEX_t), compareIndices)
212213

213214
cpdef set toSet(self):
214215
cdef:

‎drivers/testDistOp.py

+43-42
Original file line numberDiff line numberDiff line change
@@ -57,51 +57,52 @@
5757
assert not d.buildDistributedH2
5858
assert nPP.kernel.horizon.value < np.inf
5959

60-
if nPP.domain == 'gradedInterval':
61-
if d.horizonToMeshSize <= 0. or nPP.kernel.horizon.value == np.inf:
62-
h = 0.03/2**(d.noRef-6)
63-
else:
64-
h = nPP.kernel.horizon.value/d.horizonToMeshSize
65-
mesh, _ = nonlocalMeshFactory(nPP.domain,
66-
kernel=nPP.kernel,
67-
boundaryCondition=HOMOGENEOUS_DIRICHLET,
68-
h=h)
69-
elif nPP.domain == 'disc':
70-
if d.horizonToMeshSize <= 0. or nPP.kernel.horizon.value == np.inf:
71-
h = 0.04/2**(nPP.noRef-3)
72-
else:
73-
h = nPP.kernel.horizon.value/d.horizonToMeshSize/np.sqrt(2)
74-
mesh, _ = nonlocalMeshFactory(nPP.domain,
75-
kernel=nPP.kernel,
76-
boundaryCondition=HOMOGENEOUS_DIRICHLET,
77-
h=h,
78-
max_volume=h**2/2,
79-
projectNodeToOrigin=False)
80-
elif d.domain == 'gradedDisc':
81-
if d.horizonToMeshSize <= 0. or nPP.kernel.horizon.value == np.inf:
82-
h = 0.06/2**(d.noRef-6)
60+
with d.timer('set up mesh and dofmap'):
61+
if nPP.domain == 'gradedInterval':
62+
if d.horizonToMeshSize <= 0. or nPP.kernel.horizon.value == np.inf:
63+
h = 0.03/2**(d.noRef-6)
64+
else:
65+
h = nPP.kernel.horizon.value/d.horizonToMeshSize
66+
mesh, _ = nonlocalMeshFactory(nPP.domain,
67+
kernel=nPP.kernel,
68+
boundaryCondition=HOMOGENEOUS_DIRICHLET,
69+
h=h)
70+
elif nPP.domain == 'disc':
71+
if d.horizonToMeshSize <= 0. or nPP.kernel.horizon.value == np.inf:
72+
h = 0.04/2**(nPP.noRef-3)
73+
else:
74+
h = nPP.kernel.horizon.value/d.horizonToMeshSize/np.sqrt(2)
75+
mesh, _ = nonlocalMeshFactory(nPP.domain,
76+
kernel=nPP.kernel,
77+
boundaryCondition=HOMOGENEOUS_DIRICHLET,
78+
h=h,
79+
max_volume=h**2/2,
80+
projectNodeToOrigin=False)
81+
elif d.domain == 'gradedDisc':
82+
if d.horizonToMeshSize <= 0. or nPP.kernel.horizon.value == np.inf:
83+
h = 0.06/2**(d.noRef-6)
84+
else:
85+
h = nPP.kernel.horizon.value/d.horizonToMeshSize
86+
mesh, _ = nonlocalMeshFactory(nPP.domain,
87+
kernel=nPP.kernel,
88+
boundaryCondition=HOMOGENEOUS_DIRICHLET,
89+
h=h,
90+
max_volume=h**2/2)
8391
else:
84-
h = nPP.kernel.horizon.value/d.horizonToMeshSize
85-
mesh, _ = nonlocalMeshFactory(nPP.domain,
86-
kernel=nPP.kernel,
87-
boundaryCondition=HOMOGENEOUS_DIRICHLET,
88-
h=h,
89-
max_volume=h**2/2)
90-
else:
91-
if d.horizonToMeshSize <= 0. or nPP.kernel.horizon.value == np.inf:
92-
mesh = nPP.mesh
93-
for _ in range(nPP.noRef):
94-
mesh = mesh.refine()
92+
if d.horizonToMeshSize <= 0. or nPP.kernel.horizon.value == np.inf:
93+
mesh = nPP.mesh
94+
for _ in range(nPP.noRef):
95+
mesh = mesh.refine()
96+
else:
97+
mesh = nPP.mesh
98+
while d.horizonToMeshSize > np.around(nPP.kernel.horizon.value/mesh.h, 5):
99+
mesh = mesh.refine()
100+
if nPP.boundaryCondition == HOMOGENEOUS_DIRICHLET:
101+
dm = dofmapFactory(nPP.element, mesh, nPP.domainIndicator)
95102
else:
96-
mesh = nPP.mesh
97-
while d.horizonToMeshSize > np.around(nPP.kernel.horizon.value/mesh.h, 5):
98-
mesh = mesh.refine()
99-
if nPP.boundaryCondition == HOMOGENEOUS_DIRICHLET:
100-
dm = dofmapFactory(nPP.element, mesh, nPP.domainIndicator)
101-
else:
102-
dm = dofmapFactory(nPP.element, mesh, nPP.domainIndicator+nPP.fluxIndicator)
103+
dm = dofmapFactory(nPP.element, mesh, nPP.domainIndicator+nPP.fluxIndicator)
103104

104-
assert d.comm.allreduce(dm.num_dofs, op=MPI.MAX) == dm.num_dofs
105+
assert d.comm.allreduce(dm.num_dofs, op=MPI.MAX) == dm.num_dofs
105106

106107
info = d.addOutputGroup('info')
107108
info.add('Global mesh', dm.mesh)

‎fem/PyNucleus_fem/mesh.py

+1-1
Original file line numberDiff line numberDiff line change
@@ -311,7 +311,7 @@ def squareWithInteractions(ax, ay, bx, by,
311311
preserveLinesVertical=[],
312312
**kwargs):
313313
if h is None:
314-
h = horizon
314+
h = horizon-1e-8
315315
if innerRadius > 0:
316316
uniform = False
317317
if not uniform:

‎fem/PyNucleus_fem/vector_{SCALAR}.pxi

+4-1
Original file line numberDiff line numberDiff line change
@@ -459,7 +459,7 @@ cdef class {SCALAR_label_lc_}multi_fe_vector:
459459
else:
460460
v1 = self
461461
alpha = other
462-
v2 = {SCALAR_label_lc_}multi_fe_vector(np.empty((v1.data.shape[0], v1.data.shape[0]), dtype={SCALAR}), v1.dm)
462+
v2 = {SCALAR_label_lc_}multi_fe_vector(np.empty((v1.data.shape[0], v1.data.shape[1]), dtype={SCALAR}), v1.dm)
463463
{SCALAR_label_lc_}assignScaled_2d(v2.data, v1.data, alpha)
464464
return v2
465465
else:
@@ -476,6 +476,9 @@ cdef class {SCALAR_label_lc_}multi_fe_vector:
476476
{SCALAR_label_lc_}assignScaled_2d(v2.data, v1.data, alpha)
477477
return v2
478478

479+
def __rmul__(self, other):
480+
return self.__mul__(other)
481+
479482
def scale(self, {SCALAR}_t[::1] other):
480483
cdef:
481484
INDEX_t i, j

‎nl/PyNucleus_nl/clusterMethodCy.pxd

+22-1
Original file line numberDiff line numberDiff line change
@@ -22,7 +22,8 @@ from PyNucleus_base.performanceLogger cimport PLogger, FakePLogger
2222
from PyNucleus_fem.DoFMaps cimport DoFMap
2323
from PyNucleus_fem.meshCy cimport meshBase
2424
from . fractionalOrders cimport fractionalOrderBase
25-
from . kernelsCy cimport FractionalKernel
25+
from . kernelsCy cimport Kernel, FractionalKernel
26+
2627

2728
cdef class transferMatrixBuilder:
2829
cdef:
@@ -67,9 +68,13 @@ cdef class tree_node:
6768
public INDEX_t distFromRoot
6869
public indexSet _dofs
6970
public indexSet _local_dofs
71+
public indexSet secondary_dofs
7072
INDEX_t _num_dofs
7173
public indexSet _cells
7274
public REAL_t[:, ::1] box
75+
public REAL_t[:, :, ::1] boxes
76+
public REAL_t[:, ::1] coords
77+
public REAL_t[::1] hVector
7378
public REAL_t[:, ::1] transferOperator
7479
public REAL_t[:, :, ::1] value
7580
public REAL_t[::1] coefficientsUp
@@ -82,6 +87,11 @@ cdef class tree_node:
8287
public REAL_t hmin
8388
public INDEX_t interpolation_order
8489
public refinementParams refParams
90+
cdef void init(self, REAL_t[:, :, ::1] boxes, REAL_t[:, ::1] coords, REAL_t[::1] hVector)
91+
cdef void releaseData(self, BOOL_t recurse=*)
92+
cdef void releaseDoFs(self, BOOL_t recurse=*)
93+
cdef void releaseLocalDoFs(self, BOOL_t recurse=*)
94+
cdef void releaseCells(self, BOOL_t recurse=*)
8595
cdef indexSet get_dofs(self)
8696
cdef indexSet get_local_dofs(self)
8797
cdef INDEX_t get_num_dofs(self)
@@ -247,3 +257,14 @@ cdef class DistributedH2Matrix_localData(LinearOperator):
247257
REAL_t[::1] x,
248258
REAL_t[::1] y) except -1
249259
cpdef tuple convert(self)
260+
261+
262+
cdef enum admissibilityType:
263+
ADMISSIBLE = 0
264+
INADMISSIBLE = 1
265+
ZERO = 2
266+
267+
268+
cpdef admissibilityType queryAdmissibility(Kernel kernel,
269+
tree_node n1,
270+
tree_node n2)

0 commit comments

Comments
 (0)