9 use hecmw,
only: kint, kreal
20 integer(kind=kint) :: eid
21 integer(kind=kint) :: etype
22 integer(kind=kint),
pointer :: nodes(:)=>null()
23 integer(kind=kint) :: n_neighbor
24 integer(kind=kint) :: n_neighbor_max
25 integer(kind=kint),
pointer :: neighbor(:)=>null()
26 real(kind=kreal) :: reflen
27 real(kind=kreal) :: xavg(3)
28 real(kind=kreal) :: dmax
29 integer(kind=kint) :: bktid
32 integer(kind=kint),
parameter,
private :: debug = 0
39 integer(kind=kint),
intent(in) :: eid
40 integer(kind=kint),
intent(in) :: etype
41 integer(kind=kint),
intent(in) :: nsurf
43 integer(kind=kint) :: n, outtype, nodes(100)
45 call getsubface( etype, nsurf, outtype, nodes )
48 allocate( surf%nodes(n) )
49 surf%nodes(1:n)=nodes(1:n)
52 surf%n_neighbor_max = 0
62 if(
associated(surf%nodes) )
deallocate( surf%nodes )
63 if(
associated(surf%neighbor) )
deallocate( surf%neighbor )
68 integer(kind=kint),
intent(in) :: file
70 integer(kind=kint) :: i
71 write(file,*)
"Element:",surf%eid,
"Surface type=",surf%etype
72 if(
associated( surf%nodes ) ) &
73 write(file,*) ( surf%nodes(i),i=1,
size(surf%nodes) )
74 if(
associated( surf%neighbor ) ) &
75 write(file,*) ( surf%neighbor(i),i=1,surf%n_neighbor )
85 integer(kind=kint) :: i, j, ii,jj, nd1, nd2, nsurf
86 integer(kind=kint) :: k, oldsize, newsize
87 integer(kind=kint),
pointer :: dumarray(:) => null()
88 integer(kind=kint) :: bktID, ncand, js
89 integer(kind=kint),
allocatable :: indexSurf(:)
90 if (debug >= 1)
write(0,*)
'DEBUG: find_surface_neighbor: start'
100 if (ncand == 0) cycle
101 allocate(indexsurf(ncand))
106 if(
associated(surf(i)%neighbor) )
then
107 if ( any( surf(i)%neighbor(1:surf(i)%n_neighbor)==j ) ) cycle
109 do ii=1,
size(surf(i)%nodes)
110 nd1 = surf(i)%nodes(ii)
111 do jj=1,
size(surf(j)%nodes)
112 nd2 = surf(j)%nodes(jj)
114 surf(i)%n_neighbor = surf(i)%n_neighbor+1
116 if( .not.
associated(surf(i)%neighbor) )
then
119 surf(i)%neighbor(1) = j
120 else if( surf(i)%n_neighbor > surf(i)%n_neighbor_max )
then
121 oldsize = surf(i)%n_neighbor_max
122 newsize = oldsize * 2
123 dumarray => surf(i)%neighbor
124 allocate( surf(i)%neighbor(newsize) )
125 surf(i)%n_neighbor_max = newsize
127 surf(i)%neighbor(k) = dumarray(k)
129 surf(i)%neighbor(oldsize+1) = j
130 deallocate( dumarray )
132 surf(i)%neighbor(surf(i)%n_neighbor) = j
140 deallocate(indexsurf)
144 if (debug >= 1)
write(0,*)
'DEBUG: find_surface_neighbor: end'
151 real(kind=kreal),
intent(in) :: cpos(2)
153 integer(kind=kint) :: i
154 real(kind=kreal) :: maxv(3)
156 if( .not.
associated(surf%neighbor) )
return
158 select case(surf%etype)
160 if( all(cpos>0.d0) .and. all(cpos<1.d0) )
return
161 if(
size(surf%neighbor)<3 )
return
162 if( all(cpos(:)>0.d0) )
return
164 if( cpos(i)< 0.d0 ) maxv(i) = -cpos(i)
168 if( all(cpos>-1.d0) .and. all(cpos<1.d0) )
return
169 if(
size(surf%neighbor)<4 )
return
170 if( cpos(1)<-1.d0 .and. dabs(cpos(2))<1.d0 )
then
172 elseif( cpos(1)>1.d0 .and. dabs(cpos(2))<1.d0 )
then
174 elseif( dabs(cpos(1))<1.d0 .and. cpos(2)<-1.d0 )
then
176 elseif( dabs(cpos(1))<1.d0 .and. cpos(2)>1.d0 )
then
178 elseif( cpos(1)<-1.d0 .and. cpos(2)<-1.d0 )
then
179 if( cpos(1)>cpos(2) )
then
184 elseif( cpos(1)<-1.d0 .and. cpos(2)>1.d0 )
then
185 if( dabs(cpos(1))>cpos(2) )
then
190 elseif( cpos(1)>1.d0 .and. cpos(2)<-1.d0 )
then
191 if( cpos(1)>dabs(cpos(2)) )
then
196 elseif( cpos(1)>1.d0 .and. cpos(2)>1.d0 )
then
197 if( cpos(1)>cpos(2) )
then
204 stop
"type of surface element not defined"
214 real(kind=kreal),
intent(in) :: coord(:)
216 integer(kind=kint) :: nn, i, j, iss
219 nn =
size(surf(j)%nodes)
221 iss = surf(j)%nodes(i)
222 elem(1:3,i) = coord(3*iss-2:3*iss)
233 real(kind=kreal),
intent(in) :: currpos(:)
234 integer(kind=kint) :: nn, i, j, iss
238 nn =
size(surf(j)%nodes)
240 iss = surf(j)%nodes(i)
241 elem(1:3,i) = currpos(3*iss-2:3*iss)
244 xmin(i) = minval(elem(i,1:nn))
245 xmax(i) = maxval(elem(i,1:nn))
246 xsum(i) = sum(elem(i,1:nn))
248 surf(j)%xavg(:) = xsum(:) / nn
249 lx(:) = xmax(:) - xmin(:)
250 surf(j)%dmax = maxval(lx) * 0.5d0
258 real(kind=kreal),
intent(in) :: x0(3)
259 real(kind=kreal),
intent(in) :: exp_rate
260 real(kind=kreal) :: dist(3), er
261 er = max(exp_rate, 1.d0)
262 dist(:) = abs(x0(:) - surf%xavg(:))
263 if ( maxval(dist(:)) < surf%dmax * er )
then
274 type(
bucketdb),
intent(inout) :: bktDB
275 real(kind=kreal) :: x_min(3), x_max(3), d_max
276 integer(kind=kint) :: nsurf, i, j
277 if (debug >= 1)
write(0,*)
'DEBUG: update_surface_bucket_info: start'
279 if (nsurf == 0)
return
280 x_min(:) = surf(1)%xavg(:)
281 x_max(:) = surf(1)%xavg(:)
285 if (surf(i)%xavg(j) < x_min(j)) x_min(j) = surf(i)%xavg(j)
286 if (surf(i)%xavg(j) > x_max(j)) x_max(j) = surf(i)%xavg(j)
287 if (surf(i)%dmax > d_max) d_max = surf(i)%dmax
291 x_min(j) = x_min(j) - d_max
292 x_max(j) = x_max(j) + d_max
305 if (debug >= 1)
write(0,*)
'DEBUG: update_surface_bucket_info: end'
This module provides bucket-search functionality It provides definition of bucket info and its access...
subroutine, public bucketdb_allocate(bktdb)
Allocate memory before actually registering members Before allocating memory, bucketDB_registerPre ha...
subroutine, public bucketdb_registerpre(bktdb, bid)
Pre-register for just counting members to be actually registered Bucket ID has to be obtained with bu...
integer(kind=kint) function, public bucketdb_getbucketid(bktdb, x)
Get bucket ID that includes given point.
subroutine, public bucketdb_register(bktdb, bid, sid)
Register member Before actually register, bucketDB_allocate has to be called.
integer(kind=kint) function, public bucketdb_getnumcand(bktdb, bid)
Get number of candidates within neighboring buckets of a given bucket Bucket ID has to be obtained wi...
subroutine, public bucketdb_getcand(bktdb, bid, ncand, cand)
Get candidates within neighboring buckets of a given bucket Number of candidates has to be obtained w...
subroutine, public bucketdb_setup(bktdb, x_min, x_max, dmin, n_tot)
Setup basic info of buckets.
This module encapsulate the basic functions of all elements provide by this software.
integer(kind=kind(2)) function getnumberofnodes(etype)
Obtain number of nodes of the element.
subroutine getelementcenter(fetype, localcoord)
Return natural coordinate of the center of surface element.
integer, parameter fe_tri6n
integer(kind=kind(2)) function getnumberofsubface(etype)
Obtain number of sub-surface.
subroutine getsubface(intype, innumber, outtype, nodes)
Find the definition of surface of the element.
integer, parameter fe_tri3n
integer, parameter fe_quad4n
real(kind=kreal) function getreferencelength(fetype, nn, localcoord, elecoord)
This function calculates reference length at a point in surface.
integer, parameter fe_quad8n
This module provides aux functions.
This module manage surface elements in 3D It provides basic definition of surface elements (trianglar...
subroutine initialize_surf(eid, etype, nsurf, surf)
Initializer.
integer(kind=kint), parameter l_max_elem_node
subroutine update_surface_reflen(surf, coord)
Compute reference length of surface elements.
subroutine write_surf(file, surf)
Write out elemental surface.
subroutine update_surface_box_info(surf, currpos)
Update info of cubic box including surface elements.
subroutine find_surface_neighbor(surf, bktdb)
Find neighboring surface elements.
integer(kind=kint), parameter l_max_elem_surf
integer(kind=kint), parameter n_neighbor_max_init
integer(kind=kint), parameter l_max_surface_node
subroutine finalize_surf(surf)
Memeory management subroutine.
logical function is_in_surface_box(surf, x0, exp_rate)
Check if given point is within cubic box including surface element.
subroutine update_surface_bucket_info(surf, bktdb)
Update bucket info for searching surface elements.
integer(kind=kint) function next_position(surf, cpos)
Tracing next contact position.
Structure for bucket search.
Structure to define surface group.