60 type(
mesh),
intent(inout) :: m
63 &
"mesh_interfaces : define_interfaces" 78 &
" Number of interfaces =", m%nbItf
91 type(
mesh),
intent(inout) :: m
93 integer,
dimension(:,:),
allocatable :: tab, clTab
94 integer,
dimension(:) ,
allocatable :: clType
95 integer,
dimension(CELL_MAX_NBNODES) :: cl_nodes
97 integer :: cpt1, cpt2, cpt3
98 integer :: ii, cl, idx, j1, j2, nNd
101 &
"mesh_interfaces : define_domain_boundary" 104 &
"mesh_interfaces: define_domain_boundary: mesh 'm' not valid")
118 if (tab(1,ii)==0) cycle
157 if (tab(1,ii)==0) cycle
161 if (tab(3,ii)/=0) cycle
171 j1 = m%clToNd%row(cl)
172 j2 = m%clToNd%row(cl+1)
174 cl_nodes(1:nnd) = m%clToNd%col(j1: j2-1)
180 select case(m%clType(cl))
184 cltab(1,cpt3) = cl_nodes(idx)
197 call quit(
"mesh_interfaces: define_domain_boundary:& 198 & cell type not implemented: "&
213 &
" Mesh domain boundary =", cpt1+cpt2,
" cells" 244 type(
mesh) ,
intent(inout) :: m
245 integer,
dimension(:,:),
allocatable :: tab
247 integer :: itf, itf2, itf_idx, itf_nb_vtx
248 integer :: j1, j2, jj, ll
249 integer :: dim, cl, cl_type, wdt
251 integer,
dimension(CELL_MAX_FC_NBVTX) :: itf_vtx
252 integer,
dimension(:),
allocatable :: t0, t1, t2
255 &
"mesh_interfaces : boundary_interface_scan " 260 &
"mesh_interfaces: boundary_interface_scan: mesh 'm' not valid")
269 allocate(t0(wdt), t1(wdt), t2(wdt))
275 do itf=1, m%itfToCl%nl
276 j1 = m%itfToCl%row(itf)
277 j2 = m%itfToCl%row(itf+1)
282 cl = m%itfToCl%col(j1)
292 do itf=1, m%itfToCl%nl
293 if (tab(1, itf) == 0) cycle
298 cl_type = m%clType(cl)
303 j1 = m%clToItf%row(cl)
304 j2 = m%clToItf%row(cl+1) - 1
305 loop_1:
do jj= j1, j2
306 itf2 = m%clToItf%col(jj)
308 if (itf == itf2)
then 313 tab(2, itf) = itf_idx
319 j1 = m%clToNd%row(cl) - 1
321 itf_vtx(1) = m%clToNd%col(j2)
326 j1 = m%clToNd%row(cl) - 1
329 itf_vtx(1) = m%clToNd%col(j2)
332 itf_vtx(2) = m%clToNd%col(j2)
337 itf_vtx(1:itf_nb_vtx) = &
340 j1 = m%clToNd%row(cl) - 1
343 j2 = j1 + itf_vtx(jj)
344 itf_vtx(jj) = m%clToNd%col(j2)
354 call getrow(jj, t0, wdt, m%ndToCl, itf_vtx(1))
357 call getrow(j1, t1, wdt, m%ndToCl, itf_vtx(1))
358 call getrow(j2, t2, wdt, m%ndToCl, itf_vtx(2))
363 call getrow(j1, t1, wdt, m%ndToCl, itf_vtx(1))
365 loop_2:
do ll=2, itf_nb_vtx
366 call getrow(j2, t2, wdt, m%ndToCl, itf_vtx(ll))
371 if (ll==itf_nb_vtx)
exit loop_2
380 if ((jj==0) .OR. (jj>2))
call quit(&
381 &
"mesh_interfaces: boundary_interface_scan: 3")
393 j1 = m%clType(tab(3, itf) )
395 if (j2 /= dim-1)
call quit(&
396 &
"mesh_interfaces: boundary_interface_scan: 4")
407 integer,
dimension(:,:),
allocatable :: edTab
408 integer,
dimension(:) ,
allocatable :: nnz
409 integer :: cl, j1, j2
411 integer :: cl_t, nb_vtx, nb_ed, width, nb_itf
415 &
"mesh_tools : define_2DMesh_edges" 435 m%clToItf =
graph(nnz, m%nbCl, 0)
440 call allocmem(edtab, 3, m%clToItf%nnz)
441 width = m%ndToCl%width
444 if ( m%cell_count(cl_t) == 0 ) cycle
456 m%clToItf%nc = nb_itf
462 j1 = m%clToItf%row(cl)
463 j2 = m%clToItf%row(cl+1)-1
465 m%clToItf%col(j1:j2) = edtab(3,j1:j2)
472 &
" 2D-cells --> edges graph, CPU = ", &
484 &
" edges --> 2D-cells graph, CPU = ", &
496 integer,
dimension(nb_vtx) :: vtx_nbCl
497 integer,
dimension(width, nb_vtx) :: vtx_cl
498 integer,
dimension(nb_vtx) :: cl_vtx
499 integer,
dimension(3) :: coBord
500 integer :: ed, v1, v2
507 if (l1 /= cl_t ) cycle
511 j1 = m%clToNd%row(cl)
513 cl_vtx = m%clToNd%col( j1:j2 )
520 j1 = m%ndToCl%row(l1)
521 j2 = m%ndToCl%row(l1+1)
527 vtx_cl(1:l1, v1) = m%ndToCl%col(j1:j2-1)
533 j1 = m%clToItf%row(cl)
545 & vtx_cl(1:l1, v1), l1, vtx_cl(1:l2, v2), l2)
557 if (cobord(1)==cl)
then 561 edtab(2, j1) = cobord(2)
569 l1 =
cell_dim( m%clType(cobord(1)) )
572 edtab(2, j1) = cobord(1)
585 l1 =
cell_dim( m%clType(cobord(1)) )
588 if (cobord(2)==cl)
then 589 edtab(2, j1) = cobord(3)
591 edtab(2, j1) = cobord(2)
596 l1 =
cell_dim( m%clType(cobord(2)) )
599 if (cobord(1)==cl)
then 600 edtab(2, j1) = cobord(3)
602 edtab(2, j1) = cobord(1)
607 if (cobord(1)==cl)
then 608 edtab(2, j1) = cobord(2)
610 edtab(2, j1) = cobord(1)
634 integer,
intent(out) :: cpt
636 integer :: jj, cl, cl2, ll
639 do jj=1, m%clToItf%nnz
644 if ( cl <= cl2 )
then 650 ll = m%clToItf%row(cl2)
652 do while( edtab(2,ll) /= cl )
655 edtab(3,jj) = edtab(3,ll)
669 integer,
dimension(:,:),
allocatable :: fcTab
670 integer,
dimension(:) ,
allocatable :: nnz
671 integer :: cl, j1, j2
673 integer :: cl_t, nb_vtx, nb_fc, width, nb_itf
677 &
"mesh_tools : define_3DMesh_faces" 697 m%clToItf =
graph(nnz, m%nbCl, 0)
702 call allocmem(fctab, 3, m%clToItf%nnz)
703 width = m%ndToCl%width
706 if ( m%cell_count(cl_t) == 0 ) cycle
718 m%clToItf%nc = nb_itf
724 j1 = m%clToItf%row(cl)
725 j2 = m%clToItf%row(cl+1)-1
727 m%clToItf%col(j1:j2) = fctab(3,j1:j2)
734 &
" cells --> faces graph, CPU = ", &
746 &
" faces --> cells graph, CPU = ", &
758 integer,
dimension(nb_vtx) :: vtx_nbCl
759 integer,
dimension(width, nb_vtx) :: vtx_cl
760 integer,
dimension(nb_vtx) :: cl_vtx
761 integer,
dimension(width) :: coBord, aux
762 integer :: fc, v1, v2
763 integer :: l1, l2, ii
769 if (l1 /= cl_t ) cycle
773 j1 = m%clToNd%row(cl)
775 cl_vtx = m%clToNd%col( j1:j2 )
782 j1 = m%ndToCl%row(l1)
783 j2 = m%ndToCl%row(l1+1)
789 vtx_cl(1:l1, v1) = m%ndToCl%col(j1:j2-1)
795 j1 = m%clToItf%row(cl)
807 & vtx_cl(1:l1, v1), l1, vtx_cl(1:l2, v2), l2)
814 aux(1:l1) = cobord(1:l1)
817 & aux(1:l1), l1, vtx_cl(1:l2, v2), l2)
830 if (cobord(1)==cl)
then 834 fctab(2, j1) = cobord(2)
842 l1 =
cell_dim( m%clType(cobord(1)) )
845 fctab(2, j1) = cobord(1)
859 l1 =
cell_dim( m%clType(cobord(1)) )
862 if (cobord(2)==cl)
then 863 fctab(2, j1) = cobord(3)
865 fctab(2, j1) = cobord(2)
870 l1 =
cell_dim( m%clType(cobord(2)) )
873 if (cobord(1)==cl)
then 874 fctab(2, j1) = cobord(3)
876 fctab(2, j1) = cobord(1)
881 if (cobord(1)==cl)
then 882 fctab(2, j1) = cobord(2)
884 fctab(2, j1) = cobord(1)
908 integer,
intent(out) :: cpt
910 integer :: jj, cl, cl2, ll
913 do jj=1, m%clToItf%nnz
918 if ( cl <= cl2 )
then 924 ll = m%clToItf%row(cl2)
926 do while( fctab(2,ll) /= cl )
929 fctab(3,jj) = fctab(3,ll)
944 integer,
dimension(:,:),
allocatable :: vtxTab
945 integer,
dimension(:) ,
allocatable :: nnz
946 integer :: cl, j1, j2, cl_t
948 integer :: width, nb_itf
952 &
"mesh_tools : define_1DMesh_vertexes" 972 m%clToItf=
graph(nnz, m%nbCl, 0)
977 call allocmem(vtxtab, 3, m%clToItf%nnz)
978 width = m%ndToCl%width
986 m%clToItf%nc = nb_itf
992 j1 = m%clToItf%row(cl)
993 j2 = m%clToItf%row(cl+1)-1
995 m%clToItf%col(j1:j2) = vtxtab(3,j1:j2)
1002 &
" 1D-cells --> vertex graph, CPU = ", &
1014 &
" vertex --> 1D-cells graph, CPU = ", &
1026 integer,
dimension(2) :: cl_vtx
1027 integer,
dimension(width) :: coBord
1028 integer :: vtx, v, l1, l2
1036 j1 = m%clToNd%row(cl)
1038 cl_vtx = m%clToNd%col( j1:j2 )
1042 j1 = m%clToItf%row(cl)
1048 l1 = m%ndToCl%row(v)
1049 l2 = m%ndToCl%row(v+1)
1052 cobord(1:j2) = m%ndToCl%col(l1:l2-1)
1064 if (cobord(1)==cl)
then 1068 vtxtab(2, j1) = cobord(2)
1076 l1 =
cell_dim( m%clType(cobord(1)) )
1079 vtxtab(2, j1) = cobord(1)
1092 l1 =
cell_dim( m%clType(cobord(1)) )
1095 if (cobord(2)==cl)
then 1096 vtxtab(2, j1) = cobord(3)
1098 vtxtab(2, j1) = cobord(2)
1103 l1 =
cell_dim( m%clType(cobord(2)) )
1106 if (cobord(1)==cl)
then 1107 vtxtab(2, j1) = cobord(3)
1109 vtxtab(2, j1) = cobord(1)
1114 if (cobord(1)==cl)
then 1115 vtxtab(2, j1) = cobord(2)
1117 vtxtab(2, j1) = cobord(1)
1124 call quit(
"mesh_tools: define_1DMesh_vertexes: non conformal 1D mesh")
1142 integer,
intent(out) :: cpt
1144 integer :: jj, cl, cl2, ll
1147 do jj=1, m%clToItf%nnz
1152 if ( cl <= cl2 )
then 1158 ll = m%clToItf%row(cl2)
1160 do while( vtxtab(2,ll) /= cl )
1163 vtxtab(3,jj) = vtxtab(3,ll)
1186 integer,
dimension(:) ,
intent(in) :: type
1187 integer,
dimension(:,:),
intent(in) :: nodes
1189 integer,
dimension(:) ,
allocatable :: clType
1190 integer,
dimension(:,:),
allocatable :: clTab
1192 integer :: nn, j1, j2, nNd, cl
1196 &
"mesh_tools : add_cells =", nn,
" new cells" 1209 cltype(nn+1:) = m%clType
1212 j1 = m%clToNd%row(cl)
1213 j2 = m%clToNd%row(cl+1)
1216 cltab(1:nnd, cl + nn) = m%clToNd%col(j1: j2-1)
1225 cltab(:,1:nn) = nodes
1229 call allocmem(m%clType, nn + m%nbCl)
1240 if (.NOT.
valid(m))
call quit(
'mesh_interfaces: & 1241 & add_cells: construction not valid')
subroutine boundary_interface_scan(tab, m)
BUILDS Tab(1:3,1:mnbItf)
is the csr matrix well defined ?
DERIVED TYPE graph: sparse matrices of boolean in CSR format
integer, parameter cell_edg_2
Quadratic edge.
integer, parameter cell_tot_nb
Number of CELL types.
integer, parameter cell_tet
Tetrahedron.
deallocate memory for real(RP) arrays
integer, dimension(cell_tot_nb), public cell_nbvtx
Number of vertexes for each cell type.
subroutine, public mesh_clear_2(m)
Destructor for mesh type clears the data tn the mesh built after 'call create(...)'.
integer, dimension(cell_tot_nb), public cell_nbfc
Number of faces for each cell type.
subroutine add_cells(m, type, nodes)
Adds N cells to the mesh These new cells are added in first position.
subroutine, public define_interfaces(m)
Defines the mesh interfaces.
integer, parameter, public sp
simple precision for real numbers
subroutine define_3dmesh_faces(m)
integer, parameter cell_edg
Edge (line segment)
integer, dimension(cell_tot_nb), public cell_dim
Dimension for each cell type.
integer, dimension(cell_max_nbvtx, cell_max_nbfc, cell_tot_nb), public cell_fc_vtx
CELL_FC_VTX(1:n, fc, cl) = vertexes for the face fc of the cell cl, with n = CELL_FC_NBVTX(fc, cl)
REAL NUMBERS PRECISION IN CHORAL: selects simple/double/quad
integer, dimension(cell_max_nbfc, cell_tot_nb), public cell_fc_nbvtx
CELL_FC_NBVTX(fc, cl) = number of vertexes for the face fc of the cell cl.
subroutine define_1dmesh_vertexes(m)
integer, parameter cell_trg
Triangle.
subroutine number_vertexes(cpt)
subroutine number_edges(cpt)
ALGEBRAIC OPERATIONS ON SETS
integer, dimension(cell_tot_nb), public cell_nbed
Number of edges for each cell type.
DEFINE THE INTERFACES AND THE BOUNDARY CELLS OF A mesh
subroutine, public define_domain_boundary(m)
For a mesh of dimension d The boundary of the mesh domain is a collection of cells of dimension d-1...
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.
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
subroutine define_2dmesh_edges(m)
integer, parameter cell_vtx
Vertex.
subroutine, public cap_sorted_set(cpt, t, n, t1, n1, t2, n2)
t(1:cpt) = t1(1:n1) \cap t2(1:n2)
integer, dimension(2, cell_max_nbed, cell_tot_nb), public cell_ed_vtx
CELL_ED_VTX(1:2, ed, cl) = vertexes for the edge ed of the cell cl.
The type graph stores sparse matrices of boolean in CSR format.
DEFINITION OF GEOMETRICAL CELLS (for meshes)
subroutine number_faces(cpt)
integer, parameter, public cell_max_nbnodes
Cell maximal number of nodes.