31 #include "mmg/mmg2d/libmmg2df.h" 33 public :: mmg2d_retrievemesh, mmg2d_passmesh
54 type(
mesh),
intent(inout) :: m
55 real(DP) ,
dimension(:) :: val
57 integer ,
optional,
intent(in) :: verb
58 real(DP),
optional,
intent(in) :: hausd
59 real(DP),
optional,
intent(in) :: hMax
63 mmg5_data_ptr_t :: mmgmesh, mmgsol
66 &
"mesh_mmg : reMesh_levelSet" 68 if (.NOT.
valid(m))
call quit(
"mesh_mmg: reMesh_levelSet:& 69 & mesh 'm' not valid")
71 if (
size(val,1) /= m%nbNd)
call quit(
"mesh_mmg: reMesh_levelSet:& 72 & wrong size for 'val'")
79 &
"mesh_mmg: reMesh_levelSet:& 80 & dimension 2 mesh withour triangles")
87 call mmg2d_passmesh(mmgmesh, mmgsol, m)
90 call mmg2d_set_solsize(mmgmesh, mmgsol, mmg5_vertex, m%nbNd, mmg5_scalar, ier)
92 call mmg2d_set_scalarsol(mmgsol, val(ii), ii, ier)
105 if (
present(verb) )
then 106 call mmg2d_set_iparameter(mmgmesh,mmgsol,&
107 & mmg2d_iparam_verbose, verb, ier)
109 call mmg2d_set_iparameter(mmgmesh,mmgsol,&
110 & mmg2d_iparam_verbose, 0, ier)
116 if (
present(hausd) )
then 117 call mmg2d_set_dparameter(mmgmesh, mmgsol, &
118 & mmg2d_dparam_hausd, hausd, ier)
120 call mmg2d_set_dparameter(mmgmesh, mmgsol, &
121 & mmg2d_dparam_hausd, 0.01_dp, ier)
126 if (
present(hmax) )
then 127 call mmg2d_set_dparameter(mmgmesh, mmgsol, &
128 & mmg2d_dparam_hmax, hmax, ier)
133 call mmg2d_set_iparameter(mmgmesh,mmgsol,mmg2d_iparam_iso,1,ier)
139 call mmg2d_mmg2dls(mmgmesh, mmgsol, ier)
149 call mmg2d_retrievemesh(m, mmgmesh, mmgsol)
153 call mmg2d_free_all(mmg5_arg_start, &
154 mmg5_arg_ppmesh,mmgmesh,mmg5_arg_ppmet,mmgsol, &
159 if ( m%cell_Count(
cell_tet)==0 )
call quit(
"mesh_mmg: reMesh_levelSet:& 160 & dimension 3 mesh withour tetrahedra")
162 call quit(
'mesh_mmg: reMesh_levelSet:& 166 call quit(
"mesh_mmg: reMesh_levelSet: wrong dimension")
174 call quit(
"mesh_mmg: reMesh_levelSet: set 'WITH_MMG=1' in the& 175 & file 'CONFIG.in', re-configure (mmg must be installed)")
190 subroutine mmg2d_passmesh(mmgMesh, mmgSol, m1)
191 mmg5_data_ptr_t,
intent(out) :: mmgmesh, mmgsol
192 type(
mesh) ,
intent(in) :: m1
194 integer :: ier, nNd, nEd, nTr, cl, ii, j1
197 &
"mesh_mmg : mmg2D_passMesh" 202 call mmg2d_init_mesh(mmg5_arg_start, &
203 mmg5_arg_ppmesh, mmgmesh, mmg5_arg_ppmet, mmgsol, &
210 call mmg2d_set_meshsize(mmgmesh, nnd, ntr, ned, ier)
215 call mmg2d_set_vertex(mmgmesh, &
216 &m1%nd(1,ii), m1%nd(2,ii), 0, ii, ier)
225 j1 = m1%clToNd%row(cl)
226 ii = m1%clToNd%col(j1)
228 call mmg2d_set_vertex(mmgmesh, &
229 &m1%nd(1,ii), m1%nd(2,ii), m1%clTag(cl), ii, ier)
231 call mmg2d_set_requiredvertex(mmgmesh, ii, ier)
243 j1 = m1%clToNd%row(cl)
245 call mmg2d_set_edge(mmgmesh, &
246 & m1%clToNd%col(j1), m1%clToNd%col(j1+1), &
247 & m1%clTag(cl), ii, ier)
257 j1 = m1%clToNd%row(cl)
259 call mmg2d_set_triangle(mmgmesh, &
260 & m1%clToNd%col(j1), m1%clToNd%col(j1+1), m1%clToNd%col(j1+2),&
261 & m1%clTag(cl), ii, ier)
265 end subroutine mmg2d_passmesh
267 subroutine mmg2d_retrievemesh(m2, mmgMesh, mmgSol)
268 type(
mesh) ,
intent(inout) :: m2
269 mmg5_data_ptr_t,
intent(inOut) :: mmgmesh, mmgsol
271 logical,
dimension(:) ,
allocatable :: vtx_req
272 integer,
dimension(:) ,
allocatable :: vtx_ref
273 integer,
dimension(:,:),
allocatable :: clTab
275 integer :: ier, nNd, nEd, nTr, nb_vtx_cell, ii, jj, ll
276 integer :: v1, v2, v3, req, ref
280 &
"mesh_mmg : mmg2D_retrieveMesh" 284 call mmg2d_get_meshsize(mmgmesh, nnd, ntr, ned, ier)
291 call mmg2d_get_vertex(mmgmesh, x, y, vtx_ref(ii), ll, req, ier)
298 nb_vtx_cell = nb_vtx_cell + 1
300 vtx_req(ii) = .false.
305 ii = nb_vtx_cell + ned + ntr
325 if (vtx_req(ii))
then 328 m2%clTag(jj) = vtx_ref(ii)
342 call mmg2d_get_edge(mmgmesh, v1, v2, ref, ll, req, ier)
351 jj = nb_vtx_cell + ned
357 call mmg2d_get_triangle(mmgmesh, v1, v2, v3, ref, req, ier)
369 if (.NOT.
valid(m2))
call quit(
'mesh_mmg: & 370 & mmg2d_retrieveMesh: construction not valid')
372 end subroutine mmg2d_retrievemesh
is the csr matrix well defined ?
DERIVED TYPE graph: sparse matrices of boolean in CSR format
HANDLE THE COMMUNICATION WITH THE LIBRARY MMG
integer, parameter cell_tet
Tetrahedron.
deallocate memory for real(RP) arrays
integer, parameter cell_edg
Edge (line segment)
REAL NUMBERS PRECISION IN CHORAL: selects simple/double/quad
integer, parameter cell_trg
Triangle.
allocate memory for real(RP) arrays
integer choral_verb
Verbosity level.
subroutine, public mesh_create_end(m)
mesh constructor : final step
subroutine, public mesh_create_cltab(m, clTab)
mesh constructor, second step, cell types and cell nodes are known cell types are in mclType cell nod...
DEFINITION OF GLOBAL VARIABLES FOR THE LIBRARY CHORAL
integer, parameter cell_vtx
Vertex.
DEFINITION OF GEOMETRICAL CELLS (for meshes)
subroutine, public remesh_levelset(m, val, verb, hausd, hMax)