104 real(RP),
dimension(:,:),
allocatable :: nd
107 integer ,
dimension(:) ,
allocatable :: cltype
110 integer ,
dimension(:) ,
allocatable :: cltag
114 integer,
dimension(CELL_TOT_NB) :: cell_count = 0
191 call clear(m%itfToCl)
192 call clear(m%clToItf)
206 write(*,*)
"mesh_mod : print" 208 write(*,*)
" Status = valid" 210 write(*,*)
" Status = not valid" 213 write(*,*)
" Number of nodes =", m%nbNd
216 write(*,*)
" Interfaces assembled" 217 write(*,*)
" Number of interfaces =", m%nbItf
219 write(*,*)
" Interfaces not assembled" 222 write(*,*)
" Number of cells =", m%nbCl
223 write(*,*)
" Number of cells per cell type" 226 if ( m%cell_count(ii) == 0 ) cycle
230 &
" =", m%cell_count(ii)
243 if (.NOT. (
allocated(m%nd) ))
return 244 if (.NOT. (
allocated(m%clType) ))
return 245 if (.NOT. (
allocated(m%clTag) ))
return 246 if (
valid(m%clToNd) < 0 )
return 247 if (
valid(m%ndToCl) < 0 )
return 249 b = ( m%nbCl>0 ) .AND. ( m%nbNd>0 )
250 b = b .AND. ( m%nbCl == sum(m%cell_count) )
251 b = b .AND. ( m%dim>0 ) .AND. ( m%dim<4 )
252 b = b .AND. ( all( shape(m%nd) == (/3,m%nbNd/)) )
253 b = b .AND. (
size(m%clType,1) == m%nbCl )
254 b = b .AND. (
size(m%clTag,1) == m%nbCl )
255 b = b .AND. ( m%clToNd%nl == m%nbCl )
256 b = b .AND. ( m%ndToCl%nl == m%nbNd )
260 b = b .AND. (
valid(m%clToItf) > 0 )
261 b = b .AND. (
valid(m%itfToCl) > 0 )
264 b = b .AND. ( m%clToItf%nl == m%nbCl )
265 b = b .AND. ( m%itfToCl%nl == m%nbItf)
284 function readmesh(fileName, fileFormat)
result(m)
286 character(len=*),
intent(in) :: fileName, fileFormat
291 &
"mesh_mod : readMesh = "//trim(filename)
294 inquire(file=filename, exist=b)
295 if (.NOT.b)
call quit(
'mesh_mod: readMesh: uncorrect file ')
297 open(unit=10, file=filename)
298 select case(trim(fileformat))
303 call quit(
'mesh_mod: readMesh: uncorrect file format' )
309 if (.NOT.
valid(m))
call quit(
'mesh_mod: & 310 & readMesh: construction not valid')
318 integer :: ii, jj, kk, ll, nn
319 integer,
dimension(CELL_MAX_NBNODES + 8) :: tab
320 integer,
dimension(:,:),
allocatable :: clTab
327 if (nn<=0)
call quit(
"mesh_mod: gmesh_readMesh: 1" )
333 read (10, *) jj, m%nd(1:3,ii)
337 read(10,*) ;
read(10,*)
339 if (nn<=0)
call quit(
"mesh_mod: gmesh_readMesh: 2" )
348 read(10,*) jj, kk, ll, m%clTag(ii)
364 cltab(1:jj, ii) = tab(kk:ll)
381 integer ,
dimension(:,:) :: clTab
383 integer,
dimension(:),
allocatable :: nnz
384 integer :: ii, jj, ll
389 do ii=1,
size(cltab, 2)
398 m%clToNd =
graph(nnz, m%nbCl, m%nbNd)
403 do ii=1,
size(cltab, 2)
406 m%clToNd%col(jj: jj+ll-1) = cltab(1:ll, ii)
419 subroutine mesh_write(m, fileName, fileFormat, cell_tags)
421 character(len=*) ,
intent(in) :: fileName, fileFormat
422 logical,
optional,
intent(in) :: cell_tags
425 &
"mesh_mod : writee = "//trim(filename)
427 if (.NOT.
valid(m))
call quit(
'mesh_mod: mesh_write: mesh not valid')
429 open(unit=10, file=filename)
430 select case(trim(fileformat))
435 call quit(
'mesh_mod: mesh_write: uncorrect file format' )
444 logical,
optional,
intent(in) :: cell_tags
446 integer :: ii, tt, j1, j2
448 write(10,
'(A11)')
"$MeshFormat" 449 write(10,
'(A7)')
"2.2 0 8" 450 write(10,
'(A14)')
"$EndMeshFormat" 451 write(10,
'(A6)')
"$Nodes" 454 write(10, *) ii, m%nd(1:3,ii)
456 write(10,
'(A9)')
"$EndNodes" 457 write(10,
'(A9)')
"$Elements" 465 j2 = m%clToNd%row(ii+1) - 1
467 write(10, *) ii, tt, 2, m%clTag(ii), 1, m%clToNd%col(j1:j2)
471 write(10,
'(A12)')
"$EndElements" 473 if (
present(cell_tags))
then 476 write(10,fmt=
"(A)")
'$ElementData' 477 write(10,fmt=
"(A)")
'1' 478 write(10,fmt=
"(A)")
"Cell tags" 479 write(10,fmt=
"(A)")
'1' 481 write(10,fmt=
"(A)")
'3' 486 write(10,*) ii, m%clTag(ii)
488 write(10,fmt=
"(A)")
'$EndElementData' 514 real(RP),
intent(in) :: a, b
515 integer ,
intent(in) :: n
519 integer,
dimension(:),
allocatable :: p
521 if (
choral_verb>0)
write(*,*)
"mesh_mod : create_1D" 524 h = (b-a)/
real(n, RP) 531 xi = a +
real(ii,
rp)*h
532 m%nd(:,ii+1) = (/ xi, 0._rp, 0._rp/)
547 m%clToNd =
graph(p, m%nbCl, m%nbNd)
552 m%clToNd%col(j1:j1+1) = (/ii, ii+1/)
558 if (.NOT.
valid(m))
call quit(
'mesh_mod: & 559 & mesh_create_1D: construction not valid')
569 integer :: cl, ct, dim
571 if (
choral_verb>0)
write(*,*)
"mesh_mod : create_end =", &
572 &
size(m%clType,1),
" cells " 574 m%nbCl =
size(m%clType,1)
575 m%nbNd =
size(m%nd,2)
585 m%cell_count(ct) = m%cell_count(ct) + 1
588 if (dim>m%dim) m%dim = dim
598 subroutine getndcoord(ndCoord, s, ndInd, n, m)
599 integer ,
intent(in) :: n, s
600 real(RP),
dimension(3,s),
intent(out) :: ndCoord
601 integer ,
dimension(n) ,
intent(in) :: ndInd
602 type(
mesh) ,
intent(in) :: m
607 ndcoord(:,jj) = m%nd( :, ndind(jj) )
622 logical,
dimension(:),
allocatable :: flag
623 type(
mesh) ,
intent(in) :: m
624 integer ,
intent(in) :: dim
625 procedure(r3tor) ,
optional :: f
627 integer ,
dimension( CELL_MAX_NBNODES) :: p
628 real(RP),
dimension(3, CELL_MAX_NBNODES) :: X
629 integer :: cl, nbNd, ii, cpt, ct
634 &
"mesh_mod : flag_mesh_cells" 636 &
'mesh_mod: flag_mesh_cells: mesh not valid')
652 if (.NOT.flag(cl)) cycle
659 loop_nodes:
do ii=1, nbnd
664 if (.NOT.bool)
exit loop_nodes
669 if (.NOT.bool) cpt = cpt - 1
677 &
" Number of flaged cells =", cpt
698 integer ,
dimension(:) ,
allocatable :: V
699 real(RP) ,
intent(out) :: dist
700 real(RP),
dimension(:,:),
intent(in) :: pt
701 type(
mesh) ,
intent(in) :: m
702 type(
graph) ,
intent(in) :: ndToNd
703 integer ,
optional ,
intent(in) :: start
705 integer ,
dimension( ndToNd%width) :: p
706 real(RP),
dimension( ndToNd%width) :: d
707 real(RP),
dimension(3,ndToNd%width) :: Y
708 real(RP),
dimension(3) :: x
709 integer ,
dimension(1) :: ind
711 integer :: ii, jj, nd, cpt, wdt, sz
714 if (
size(pt,1)/=3)
call quit(
"mesh_mod:& 715 & closest_node: wrong argument 'pt'" )
724 if (
present(start)) nd = start
730 a = sum( (m%nd(:,nd) - x)**2)
737 call getrow(sz, p, wdt, ndtond, nd)
743 y(:,jj) = y(:,jj)-x(:)
745 y(:,1:sz) = y(:,1:sz)**2
747 d(jj) = sum( y(:,jj) )
749 ind = minloc( d(1:sz) )
767 a =
re(cpt,
size(pt,2))
769 &
"mesh_mod : closest_node = & 770 &averagge neighbour nodes used =",
real(a,4)
778 type(
mesh),
intent(in) :: m
779 integer ,
intent(in) :: cl
781 integer,
dimension(:),
allocatable :: p
783 real(RP),
dimension(3) :: AB, AC, v
788 select case(m%clType(cl))
793 call getrow(sz, p, wdt, m%clToNd, cl)
794 ab = m%nd(:,p(2)) - m%nd(:,p(1))
801 call getrow(sz, p, wdt, m%clToNd, cl)
802 ab = m%nd(:,p(2)) - m%nd(:,p(1))
803 ac = m%nd(:,p(3)) - m%nd(:,p(1))
810 call quit(
'mesh_mod: cell_orientation:& 811 & cell type not supported' )
824 integer ,
intent(out) :: nb_itf
825 integer ,
dimension(n),
intent(out) :: eps
826 type(
mesh),
intent(in) :: m
827 integer ,
intent(in) :: cl, n
829 integer,
dimension(2) :: p2
830 integer,
dimension(CELL_MAX_NBITF) :: p1
832 integer :: ii, coBound
835 if (m%nbItf<=0)
call quit( &
836 &
"mesh_mod: interface_orientation:& 837 & try 'call define_interfaces(m)'" )
848 call getrow(cobound, p2, 2, m%itfToCl, p1(ii))
852 if (p2(1)==cl) eps(ii) = 1
type(mesh) function readmesh(fileName, fileFormat)
Constructor for the type mesh
DEFINITION OF ABSTRACT_INTERFACES FOR THE LIBRARY CHORAL
subroutine mesh_clear(m)
Destructor for mesh type.
integer, dimension(:), allocatable, public gmsh_to_cell
Dictionnary GMSH cell desc. –> CHORAL cell desc.
DERIVED TYPE graph: sparse matrices of boolean in CSR format
integer, parameter cell_edg_2
Quadratic edge.
subroutine mesh_print(m)
display cell types in the mesh
integer, parameter cell_tot_nb
Number of CELL types.
deallocate memory for real(RP) arrays
subroutine, public mesh_clear_2(m)
Destructor for mesh type clears the data tn the mesh built after 'call create(...)'.
subroutine gmesh_writemesh(m, cell_tags)
write mesh to gmsh ASCII file
subroutine mesh_write(m, fileName, fileFormat, cell_tags)
Write mesh 'm' to file 'fileName' with format 'fileFormat'.
integer, parameter, public rp
real(kind=RP) = real precision in the code REAL_TOL = epsilon to test real equality ...
subroutine, public flag_mesh_cells(flag, m, dim, f)
OUTPUT: flag(:) array of logical with size mnbCl.
integer, parameter cell_edg
Edge (line segment)
integer, dimension(cell_tot_nb), public cell_dim
Dimension for each cell type.
subroutine, public interface_orientation(nb_itf, eps, n, m, cl)
eps(ii) = +1 if the interface ii of the cell cl in mesh m is oriented outwards the cell cl = -1 other...
conversion integers or rational to real
REAL NUMBERS PRECISION IN CHORAL: selects simple/double/quad
logical function mesh_valid(m)
Verify allocation.
integer, dimension(cell_tot_nb), public cell_nbnodes
Number of nodes for each cell type.
integer, parameter cell_trg
Triangle.
subroutine gmesh_readmesh(m)
read mesh from gmsh ASCII file
integer, parameter cell_trg_2
Quadratic triangle.
integer, parameter, public cell_max_nbitf
Cell maximal number of interfaces.
subroutine, public getndcoord(ndCoord, s, ndInd, n, m)
returns node coordinates of node indexes in ndInd
allocate memory for real(RP) arrays
character(len=4), dimension(cell_tot_nb), public cell_name
CELL_ARRAY Arrays describing cells.
integer choral_verb
Verbosity level.
integer, dimension(cell_tot_nb), public cell_to_gmsh
Dictionnary cell desc. –> GMSH cell desc.
subroutine, public mesh_create_end(m)
mesh constructor : final step
real(rp), parameter, public real_tol
subroutine, public closest_node(V, dist, pt, m, nDToNd, start)
For each point X=pt(:,ii), determines a node V in the mesh m that minimises |X-V|.
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.
real(rp) function, dimension(3), public crossprod_3d(u, v)
Cross product in R**3.
type(mesh) function mesh_create_1d(a, b, n)
Constructor for the type mesh
logical function, public cell_orientation(m, cl)
has cell 'cl' the direct orientation ?
The type graph stores sparse matrices of boolean in CSR format.
DEFINITION OF GEOMETRICAL CELLS (for meshes)
integer, parameter, public cell_max_nbnodes
Cell maximal number of nodes.