The TELEMAC-MASCARET system  trunk
almesh.f
Go to the documentation of this file.
1 ! *****************
2  SUBROUTINE almesh
3 ! *****************
4 !
5  &(mesh,nom,ielm,spheri,cfg,fformat,nfic,equa,refine,nplan,npmax,
6  & nptfrx,nelmax,projection,lati0,longi0,convergence,rlevel)
7 !
8 !***********************************************************************
9 ! BIEF V7P3
10 !***********************************************************************
11 !
12 !brief ALLOCATES A BIEF_MESH MESH STRUCTURE.
13 !
14 !history J-M HERVOUET (LNHE)
15 !+ 05/02/2010
16 !+ V6P0
17 !+ EDGE-BASED STUCTURES ALWAYS ALLOCATED
18 !
19 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
20 !+ 13/07/2010
21 !+ V6P0
22 !+ Translation of French comments within the FORTRAN sources into
23 !+ English comments
24 !
25 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
26 !+ 21/08/2010
27 !+ V6P0
28 !+ Creation of DOXYGEN tags for automated documentation and
29 !+ cross-referencing of the FORTRAN sources
30 !
31 !history J-M HERVOUET (LNHE)
32 !+ 28/11/2011
33 !+ V6P2
34 !+ Calls to eleb3d and eleb3dt changed. Calls of READGEO3 and CPIKLE2
35 !+ and CPIKLE3 swapped (now KNOLG used in CPIKLE3).
36 !
37 !history J-M HERVOUET (LNHE)
38 !+ 20/07/2012
39 !+ V6P2
40 !+ Finding the original number of nodes in parallel, and completing
41 !+ KNOLG for upper planes in 3D.
42 !
43 !history J-M HERVOUET (LNHE)
44 !+ 07/02/2013
45 !+ V6P3
46 !+ Argument NELMAX removed in call to SEGBOR.
47 !
48 !history J-M HERVOUET (EDF R&D, LNHE)
49 !+ 11/03/2013
50 !+ V6P3
51 !+ Dimension of LIMVOI now set to (11,2).
52 !history R.NHEILI (Univerte de Perpignan, DALI)
53 !+ 24/02/2016
54 !+ V7P3
55 ! ALLOCATE BUFER FOR ERRORS
56 !
57 !history J-M HERVOUET (EDF R&D, LNHE)
58 !+ 15/05/2015
59 !+ V7P1
60 !+ Call to a new subroutine that checks the mesh (coordinates...).
61 !
62 !history Y AUDOUIN (LNHE)
63 !+ 25/05/2015
64 !+ V7P0
65 !+ Modification to comply with the hermes module
66 !
67 !history J-M HERVOUET (EDF R&D, LNHE)
68 !+ 30/07/2015
69 !+ V7P1
70 !+ Introducing MESH%IFAC, to replace MESH%FAC.
71 !
72 !history J,RIEHME (ADJOINTWARE)
73 !+ November 2016
74 !+ V7P2
75 !+ Replaced EXTERNAL statements to parallel functions / subroutines
76 !+ by the INTERFACE_PARALLEL
77 !
78 !history J-M HERVOUET (EDF R&D, LNHE)
79 !+ 08/09/2017
80 !+ V7P3
81 !+ New optional argument RLEVEL added. Parts of the subroutine not
82 !+ executed when RLEVEL > 0 and NPOIN, NPTFR, NELEM retrieved in
83 !+ structure MESH.
84 !
85 !history J-M HERVOUET (EDF R&D, LNHE)
86 !+ 23/09/2017
87 !+ V7P3
88 !+ In the convergence procedure, exact formulas are now used for
89 !+ dimensioning the arrays. Moreover IKLES and IPOBO, and a number
90 !+ of other arrays need not be overdimensioned.
91 !
92 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
93 !| CONVERGENCE |-->| IF YES, WE ARE IN A REFINEMENT PROCEDURE
94 !| EQUA |-->| NAME IN 20 CHARACTERS TO ENABLE DIFFERENT
95 !| | | OPTIONS. OPTIONS ARE:
96 !| | | "SAINT-VENANT EF"
97 !| | | "SAINT-VENANT VF"
98 !| | | "BOUSSINESQ"
99 !| IELM |-->| ELEMENT TYPE WITH THE LARGET NUMBER OF DEGREES
100 !| | | OF FREEDOM THAT WILL BE USED
101 !| LATI0 |-->| LATITUDE OF ORIGIN POINT
102 !| LONGI0 |-->| LONGITUDE OF ORIGIN POINT
103 !| MESH |-->| MESH STRUCTURE TO BE ALLOCATED
104 !| NELMAX |-->| MAXIMUM NUMBER OF ELEMENTS IN THE MESH
105 !| NFIC |-->| LOGICAL UNIT WHERE TO READ THE MESH
106 !| NOM |-->| NAME OF THE MESH
107 !| NPLAN |-->| NUMBER OF PLANES (OPTIONAL,3D MESHES OF PRISMS)
108 !| NPMAX |-->| MAXIMUM NUMBER OF POINTS IN THE MESH
109 !| NPTFRX |-->| MAXIMUM NUMBER OF BOUNDARY NODES
110 !| PROJECTION |<->| SPATIAL PROJECTION TYPE
111 !| REFINE |-->| NUMBER OF REFINEMENT LEVELS FOR CONVERGENCE
112 !| RLEVEL |-->| REFINEMENT LEVEL FOR CONVERGENCE
113 !| SPHERI |-->| LOGICAL, IF YES : SPHERICAL COORDINATES
114 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
115 !
116  USE bief, ex_almesh => almesh
118  USE interface_hermes
119 !
121  USE interface_parallel, ONLY : p_max
122  IMPLICIT NONE
123 !
124 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
125 !
126  TYPE(bief_mesh) , INTENT(INOUT) :: MESH
127  INTEGER , INTENT(IN) :: IELM
128  INTEGER , INTENT(IN) :: NFIC
129  CHARACTER(LEN=8) , INTENT(IN) :: FFORMAT
130  LOGICAL , INTENT(IN) :: SPHERI
131  CHARACTER(LEN=6) , INTENT(IN) :: NOM
132  CHARACTER(LEN=20), INTENT(IN) :: EQUA
133  INTEGER , INTENT(INOUT) :: CFG(2)
134  INTEGER , INTENT(IN) :: REFINE
135  INTEGER , INTENT(IN), OPTIONAL :: NPLAN
136  INTEGER , INTENT(IN), OPTIONAL :: NPMAX
137  INTEGER , INTENT(IN), OPTIONAL :: NPTFRX
138  INTEGER , INTENT(IN), OPTIONAL :: NELMAX,RLEVEL
139  INTEGER , INTENT(IN), OPTIONAL :: PROJECTION
140  DOUBLE PRECISION , INTENT(IN), OPTIONAL :: LATI0,LONGI0
141  LOGICAL , INTENT(IN), OPTIONAL :: CONVERGENCE
142 !
143 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
144 !
145  INTEGER D,IELM0,IELM1,STOCFG,IELB0,IELB1,NSEG,NNPMAX
146  INTEGER NNPTFRX,NNELEB,ERR,NNELMAX,NNPLAN,NPOIN,NPTFR,NELEM
147  INTEGER MXPTVS,MXELVS,NDP,IELEM,NSEGBOR,IPLAN,NPOIN_GLOB
148  INTEGER NVAR,TYP_ELEM,TYP_BND_ELEM,IERR,REFLEVEL
149  LOGICAL CONV
150 !
151 ! TEMPORARY CONNECTIVITY TABLE
152 !
153  INTEGER, ALLOCATABLE :: IKLES(:)
154 !
155 ! TEMPORARY TABLE TO NUMBER THE BOUNDARY NODES
156 !
157  INTEGER, ALLOCATABLE :: IPOBO(:)
158 !
159  INTEGER IELB0V,IELB1V,I
160 !
161  CHARACTER(LEN=80) TITLE
162 !
163  INTEGER PROJEC
164  DOUBLE PRECISION LATI,LONGI
165 !
166 ! EARTH RADIUS
167  DOUBLE PRECISION, PARAMETER :: R=6.37d6
168 !
169 !-----------------------------------------------------------------------
170 !
171  IF(PRESENT(convergence)) THEN
172  conv = convergence
173  ELSE
174  conv = .false.
175  ENDIF
176  IF(PRESENT(rlevel).AND.conv) THEN
177  reflevel=rlevel
178  ELSE
179  reflevel=0
180  ENDIF
181  IF(PRESENT(projection)) THEN
182  projec=projection
183  ELSE
184  projec=1
185  ENDIF
186 !
187 !-----------------------------------------------------------------------
188 !
189  IF(PRESENT(lati0)) THEN
190  lati=lati0
191  ELSE
192  lati=0.d0
193  ENDIF
194 !
195 !-----------------------------------------------------------------------
196 !
197  IF(PRESENT(longi0)) THEN
198  longi=longi0
199  ELSE
200  longi=0.d0
201  ENDIF
202 !
203 !-----------------------------------------------------------------------
204 !
205 ! FIRST READS THE GEOMETRY FILE TO GET NPOIN,..
206 !
207  mesh%NAME = nom
208 !
209 ! IN PARALLEL MODE, THIS IS WHERE NPTIR IS READ
210  ALLOCATE(mesh%X_ORIG)
211  ALLOCATE(mesh%Y_ORIG)
212 !
213  CALL read_mesh_info(fformat,nfic,title,nvar,npoin,typ_elem,
214  & nelem,nptfr,nptir,ndp,nnplan,
215  & mesh%X_ORIG,mesh%Y_ORIG,
216  & typ_bnd_elem,nneleb)
217 !
218 ! ALLOCATES THE TEMPORARY CONNECTIVITY TABLE
219 !
220  ALLOCATE(ikles(nelem*ndp),stat=err)
221  CALL check_allocate(err,'ALMESH:IKLES')
222 !
223 ! ALLOCATES THE TEMPORARY TABLE FOR THE BOUNDARY NODES
224 !
225  ALLOCATE(ipobo(npoin),stat=err)
226  CALL check_allocate(err,'ALMESH:IPOBO')
227 !
228 ! GETS CONNECTIVITY TABLE IKLE AND READS IPOBO IF NOT IN PARALLEL
229 !
230  CALL read_mesh_conn(fformat,nfic,npoin,typ_elem,nelem,ndp,
231  & typ_bnd_elem,nneleb,ikles,ipobo)
232 !
233 ! CALCULATES THE MAXIMUM NUMBER OF ELEMENTS AROUND A NODE MXELVS
234 ! AND THE MAXIMUM NUMBER OF SURROUNDING NODES, MXPTVS
235 !
236  CALL mxptel(mxptvs,mxelvs,ikles,ielm,
237  & npoin,nelem,ndp,ipobo,.true.)
238 !
239  DEALLOCATE(ipobo)
240 !
241 !-----------------------------------------------------------------------
242 !
243  IF(ncsize.GT.1) THEN
244 ! IN PARALLEL MODE NSEGBOR (COMPUTED IN SEGBOR) IS NOT NPTFR
245 ! IT INCLUDES INTERNAL SEGMENTS
246  CALL segbor(nsegbor,ikles,nelem,npoin)
247  ELSE
248  nsegbor=nptfr
249  ENDIF
250 !
251 !-----------------------------------------------------------------------
252 !
253 ! INITIALISES COMMONS DIMS AND NODES
254 !
255  IF(PRESENT(npmax)) THEN
256  nnpmax = npmax
257  ELSEIF(conv) THEN
258 ! EXACT FORMULA FOR ELEMENTS SUCCESSIVELY SPLIT INTO 4
259 ! UP TO LEVEL REFINE
260  nnpmax = npoin+(nelem*(4**refine-1)+nsegbor*(2**refine-1))/2
261  ELSE
262  nnpmax = npoin
263  ENDIF
264  IF(PRESENT(nptfrx)) THEN
265  nnptfrx = nptfrx
266  ELSEIF(conv) THEN
267 ! EXACT FORMULA
268  nnptfrx = nsegbor*2**refine
269  ELSE
270  nnptfrx = nptfr
271  ENDIF
272  IF(PRESENT(nelmax)) THEN
273  nnelmax = nelmax
274  ELSEIF(conv) THEN
275 ! EXACT FORMULA
276  nnelmax = nelem*4**refine
277  ELSE
278  nnelmax = nelem
279  ENDIF
280  IF(PRESENT(nplan)) THEN
281  nnplan = nplan
282  ELSE
283  nnplan = 1
284  ENDIF
285 !
286 !-----------------------------------------------------------------------
287 !
288 ! IF IN A REFINEMENT PROCEDURE, THE DIMENSIONS OF THE INITIAL MESH
289 ! ARE NOT THE CURRENT DIMENSIONS, THEY ARE ALREADY IN THE MESH
290 ! STRUCTURE.
291 !
292  IF(reflevel.GT.0) THEN
293  npoin=mesh%NPOIN
294  nptfr=mesh%NPTFR
295  nelem=mesh%NELEM
296  IF(ncsize.GT.1) THEN
297  nsegbor=nsegbor*2**reflevel
298  ELSE
299  nsegbor=nptfr
300  ENDIF
301  ENDIF
302 !
303 ! IN CALL ININDS, ALL VALUES ARE 2D VALUES
304 ! ONLY NNPLAN TELLS IF IT'S 2D OR 3D
305 !
306 ! 3D TETRAHEDRONS MESH:
307 ! ADD THE OPTIONAL ARGUMENT NNELEB
308 !
309  IF(reflevel.EQ.0) ALLOCATE(mesh%NDS(0:81,7))
310  CALL bief_ininds(npoin,nptfr,nelem,nnpmax,nnptfrx,nnelmax,
311  & nnplan,nsegbor,mesh%NDS,nneleb)
312 !
313 ! P0 AND P1 ELEMENTS
314 !
315  ielm0 = 10*(ielm/10)
316  ielm1 = ielm0 + 1
317  d = dimens(ielm0)
318 !
319 ! BOUNDARY ELEMENTS (AT THE SURFACE AND AT THE BOTTOM FOR PRISMS)
320 !
321  ielb0 = ielbor(ielm0,1)
322  ielb1 = ielbor(ielm1,1)
323 !
324 ! LATERAL BOUNDARY ELEMENTS (DIFFERENT FOR PRISMS)
325 !
326  ielb0v = ielbor(ielm0,2)
327  ielb1v = ielbor(ielm1,2)
328 !
329 !-----------------------------------------------------------------------
330 !
331 ! ALLOCATES THE ARRAYS OF REALS
332 !
333 ! COORDINATES BY ELEMENTS: XEL, YEL, ZEL
334 !
335  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%XEL)
336  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%YEL)
337  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%ZEL)
338 !
339  CALL bief_allvec(1,mesh%XEL,'XEL ',
340  & ielm0,bief_nbpel(ielm1,mesh),1,mesh)
341  CALL bief_allvec(1,mesh%YEL,'YEL ',
342  & ielm0,bief_nbpel(ielm1,mesh),1,mesh)
343 !
344  IF(d.GE.3) THEN
345  CALL bief_allvec(1,mesh%ZEL,'ZEL ',
346  & ielm0,bief_nbpel(ielm1,mesh),1,mesh)
347  ELSE
348  CALL bief_allvec(1,mesh%ZEL,'ZEL ', 0, 1,0,mesh)
349  ENDIF
350 !
351 ! SURFACES OF THE ELEMENTS: SURFAC
352 ! JAJ CAN BE USED FOR ELEMENT VOLUMES...
353 !
354  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%SURFAC)
355  CALL bief_allvec(1,mesh%SURFAC,'SURFAC',ielm0,1,1,mesh)
356 !
357 ! 1/DET : SURDET ! NOT USED IN 3D, WHY?
358 !
359  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%SURDET)
360  CALL bief_allvec(1,mesh%SURDET,'SURDET',ielm0,1,1,mesh)
361 !
362 ! LENGTHS OF THE SEGMENTS: LGSEG
363 ! CAN BE USED (IN THEORY) FOR LATERAL SURFACES IN 3D,
364 ! BUT THEN IELB0V INSTEAD OF IELB0! (2D CASE NOT AFFECTED)
365 !
366  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%LGSEG)
367  CALL bief_allvec(1,mesh%LGSEG,'LGSEG ',ielb0v,1,1,mesh)
368 !
369 ! NORMALS TO THE SEGMENTS: XSGBOR, YSGBOR, ZSGBOR
370 ! CAN BE (IN THEORY) USED FOR "NON-SIGMA" MESH FOR LATERAL NORMAL VECTORS
371 ! PER LATERAL BOUNDARY ELEMENT, BUT THEN IELB0V INSTEAD OF IELB0!
372 ! 2D CASE NOT AFFECTED
373 !
374  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%XSGBOR)
375  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%YSGBOR)
376  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%ZSGBOR)
377 ! SEE NORMAB FOR MEANING OF 4 DIMENSIONS
378  CALL bief_allvec(1,mesh%XSGBOR,'XSGBOR',ielb0v,4,1,mesh)
379  CALL bief_allvec(1,mesh%YSGBOR,'YSGBOR',ielb0v,4,1,mesh)
380  IF(d.GE.3) THEN
381  CALL bief_allvec(1,mesh%ZSGBOR,'ZSGBOR',ielb0v,4,1,mesh)
382  ELSE
383  CALL bief_allvec(1,mesh%ZSGBOR,'ZSGBOR', 0,4,0,mesh)
384  ENDIF
385 !
386 ! NORMALS AT THE NODES: XNEBOR, YNEBOR, ZNEBOR
387 !
388 ! IN 3D THEY ARE NORMAL VECTORS AT THE BOTTOM
389 ! SO THAT IELB1 REMAINS
390 !
391  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%XNEBOR)
392  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%YNEBOR)
393  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%ZNEBOR)
394 !
395  CALL bief_allvec(1,mesh%XNEBOR,'XNEBOR',ielb1,2,1,mesh)
396  CALL bief_allvec(1,mesh%YNEBOR,'YNEBOR',ielb1,2,1,mesh)
397 !
398  IF(d.GE.3) THEN !JAJ NOT USED, ACTUALLY
399  CALL bief_allvec(1,mesh%ZNEBOR,'ZNEBOR',ielb1,2,1,mesh)
400  ELSE
401  CALL bief_allvec(1,mesh%ZNEBOR,'ZNEBOR', 0,2,0,mesh)
402  ENDIF
403 !
404 ! COORDINATES BY POINTS: X, Y AND Z
405 !
406  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%X)
407  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%Y)
408  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%Z)
409  CALL bief_allvec(1,mesh%X,'X ',ielm1,1,1,mesh)
410  CALL bief_allvec(1,mesh%Y,'Y ',ielm1,1,1,mesh)
411  IF(d.GE.3) THEN
412  CALL bief_allvec(1,mesh%Z,'Z ',ielm1,1,1,mesh)
413  ELSE
414  CALL bief_allvec(1,mesh%Z,'Z ', 0,1,0,mesh)
415  ENDIF
416 !
417 ! COS AND SIN OF THE LATITUDE
418 ! WITH IELM (EXAMPLE : VELOCITY IN CORIOLIS)
419 ! COSLAT AND SINLAT ARE WORKING ARRAYS, TO WHICH
420 ! THE STRUCTURE OF X IS GIVEN TO BEGIN WITH.
421 ! THEY CAN BE EXTENDED TO ELEMENT IELM AT A LATER DATE.
422 !
423  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%COSLAT)
424  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%SINLAT)
425  IF(spheri) THEN
426  CALL bief_allvec(1,mesh%COSLAT,'COSLAT',ielm,1,2,mesh)
427  CALL bief_allvec(1,mesh%SINLAT,'SINLAT',ielm,1,2,mesh)
428  CALL cpstvc(mesh%X,mesh%COSLAT)
429  CALL cpstvc(mesh%X,mesh%SINLAT)
430  ELSE
431  CALL bief_allvec(1,mesh%COSLAT,'COSLAT', 0,1,0,mesh)
432  CALL bief_allvec(1,mesh%SINLAT,'SINLAT', 0,1,0,mesh)
433  ENDIF
434 !
435 ! DISTANCES TO BOUNDARIES : DISBOR
436 !
437  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%DISBOR)
438  CALL bief_allvec(1,mesh%DISBOR,'DISBOR',ielb0,1,1,mesh)
439 !
440 ! WORKING MATRIX (INTERNAL TO BIEF), WITH CLASSICAL STORAGE
441 !
442  stocfg = cfg(1)
443  cfg(1) = 1
444  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%M)
445  CALL bief_allmat(mesh%M,'M ',ielm,ielm,cfg,'Q','Q',mesh)
446  cfg(1) = stocfg
447 !
448 ! WORKING MATRIX BY SEGMENT
449 !
450  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%MSEG)
451 ! FROM 5.9 ON, ALWAYS DONE IN 2D (IELM=11,12,13 OR 14)
452  IF(cfg(1).EQ.3.OR.10*(ielm/10).EQ.10) THEN
453  CALL bief_allmat(mesh%MSEG,'MSEG ',ielm,ielm,cfg,'Q','Q',mesh)
454  ELSE
455  CALL bief_allmat(mesh%MSEG,'MSEG ',ielm,ielm,cfg,'0','0',mesh)
456  ENDIF
457 !
458 ! WORKING ARRAY FOR A NOT ASSEMBLED VECTOR
459 !
460  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%W)
461  CALL bief_allvec(1,mesh%W,'W ',
462  & ielm0,bief_nbpel(ielm,mesh),2,mesh)
463 !
464 ! WORKING ARRAY FOR A NORMAL VECTOR
465 !
466  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%T)
467  CALL bief_allvec(1,mesh%T,'T ',ielm,1,2,mesh)
468 !
469 ! VNOIN: ARRAY WITH NORMALS VNOIN FOR FINITE VOLUMES
470 ! CMI: COORDINATES OF THE MIDDLE OF THE SEGMENTS (KINETIC SCHEMES)
471 ! DTHAUT:
472 ! DPX,DPY: GRADIENTS OF THE BASE FUNCTIONS
473 !
474  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%VNOIN)
475  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%CMI)
476  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%AIRST)
477  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%DTHAUT)
478  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%DPX)
479  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%DPY)
480  IF(equa(1:15).EQ.'SAINT-VENANT VF') THEN
481  nseg=bief_nbseg(ielm1,mesh)
482  CALL bief_allvec(1,mesh%VNOIN ,'VNOIN ',3*nseg,1,0,mesh)
483  CALL bief_allvec(1,mesh%CMI ,'CMI ',2*nseg,1,0,mesh)
484  CALL bief_allvec(1,mesh%AIRST ,'AIRST ',2*nseg,1,0,mesh)
485  CALL bief_allvec(1,mesh%DTHAUT,'DTHAUT',ielm1 ,1,2,mesh)
486  CALL bief_allvec(1,mesh%DPX ,'DPX ',ielm0 ,3,2,mesh)
487  CALL bief_allvec(1,mesh%DPY ,'DPY ',ielm0 ,3,2,mesh)
488  ELSE
489  CALL bief_allvec(1,mesh%VNOIN ,'VNOIN ', 0,1,0,mesh)
490  CALL bief_allvec(1,mesh%CMI ,'CMI ', 0,1,0,mesh)
491  CALL bief_allvec(1,mesh%AIRST ,'AIRST ', 0,1,0,mesh)
492  CALL bief_allvec(1,mesh%DTHAUT,'DTHAUT', 0,1,0,mesh)
493  CALL bief_allvec(1,mesh%DPX ,'DPX ', 0,1,0,mesh)
494  CALL bief_allvec(1,mesh%DPY ,'DPY ', 0,1,0,mesh)
495  ENDIF
496 ! FOR COORDINATES OF CENTER OF GRAVITY OF ELEMENTS NEIGHBORING EDGES
497 !
498  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%COORDG)
499  nseg=bief_nbseg(ielm1,mesh)
500  CALL bief_allvec(1,mesh%COORDG ,'COORDG',4*nseg,1,0,mesh)
501 !
502 ! FOR PARALLEL MODE
503 !
504  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%XSEG)
505  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%YSEG)
506  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%IFAC)
507 ! THEIR ALLVEC IS IN PARINI
508  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%BUF_SEND)
509  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%BUF_RECV)
510 !
511  IF (modass .EQ.3) THEN
512  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%BUF_SEND_ERR)
513  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%BUF_RECV_ERR)
514  ENDIF
515 !
516  IF(ncsize.GT.1) THEN
517 ! XSEG
518  CALL bief_allvec(1,mesh%XSEG,'XSEG ',ielbor(ielm1,2),1,2,mesh)
519 ! YSEG
520  CALL bief_allvec(1,mesh%YSEG,'YSEG ',ielbor(ielm1,2),1,2,mesh)
521 ! IFAC
522  CALL bief_allvec(2,mesh%IFAC,'IFAC ',ielm ,1,2,mesh)
523  ELSE
524  CALL bief_allvec(1,mesh%XSEG ,'XSEG ',0,1,0,mesh)
525  CALL bief_allvec(1,mesh%YSEG ,'YSEG ',0,1,0,mesh)
526  CALL bief_allvec(2,mesh%IFAC ,'IFAC ',0,1,0,mesh)
527  ENDIF
528 !
529 !-----------------------------------------------------------------------
530 !
531 ! 1) INTEGER VALUES (ALLOCATE BECAUSE THEY ARE POINTERS)
532 !
533  IF(reflevel.EQ.0) ALLOCATE(mesh%NELEM)
534  mesh%NELEM = bief_nbpts(ielm0,mesh)
535  IF(reflevel.EQ.0) ALLOCATE(mesh%NELMAX)
536 ! WILL GIVE BACK NNELMAX FEEDED TO BIEF_ININDS ABOVE...
537  mesh%NELMAX = bief_nbmpts(ielm0,mesh)
538 !
539 !
540 ! I DO USE MESH%NPTFR FOR THE NUMBER OF LATERAL BOUNDARY NODES
541 ! IELBOR(IELM0,1) CHANGED TO IELBOR(IELM1,1)
542 ! THE PROBLEM IS, THAT FOR 3D (IELM=41):
543 ! BIEF_NBPTS(IELBOR(IELM0,1),MESH) IS THE NUMBER OF HORIZONTAL BOUNDARY ELEMENTS
544 ! BIEF_NBPTS(IELBOR(IELM0,2),MESH) IS THE NUMBER OF VERTICAL BOUNDARY ELEMENTS
545 ! BIEF_NBPTS(IELBOR(IELM1,1),MESH) IS THE NUMBER OF HORIZONTAL BOUNDARY NODES
546 ! BIEF_NBPTS(IELBOR(IELM1,2),MESH) IS THE NUMBER OF VERTICAL BOUNDARY NODES
547 !
548 ! FUNNY, BUT THE 2D CASE IS NOT AFFECTED, BECAUSE THE NUMBER OF BOUNDARY
549 ! SEGMENTS IS EQUAL TO THE NUMBER OF BOUNDARY NODES.
550 !
551  IF(reflevel.EQ.0) ALLOCATE(mesh%NPTFR)
552  mesh%NPTFR = bief_nbpts(ielbor(ielm1,2),mesh)
553  IF(reflevel.EQ.0) ALLOCATE(mesh%NPTFRX)
554  mesh%NPTFRX = bief_nbmpts(ielbor(ielm1,2),mesh)
555 !
556 ! NUMBER OF LATERAL BOUNDARY ELEMENTS
557 !
558  IF(reflevel.EQ.0) ALLOCATE(mesh%NELEB)
559  IF(reflevel.EQ.0) ALLOCATE(mesh%NELEBX)
560 !
561 ! 3D MESH
562 !
563  IF(ielm.EQ.31) THEN
564 ! HERE NNELEB HAS BEEN READ IN THE GEOMETRY FILE
565  mesh%NELEB = nneleb
566  mesh%NELEBX = nneleb
567  ELSEIF(ielm.EQ.11.OR.ielm.EQ.12.OR.ielm.EQ.13.OR.ielm.EQ.14) THEN
568  mesh%NELEB = nptfr
569  mesh%NELEBX = nnptfrx
570  ELSEIF(ielm.EQ.41) THEN
571  mesh%NELEB = nptfr*(nnplan-1)
572  mesh%NELEBX = nnptfrx*(nnplan-1)
573  ELSEIF(ielm.EQ.51) THEN
574  mesh%NELEB = 2*nptfr*(nnplan-1)
575  mesh%NELEBX = 2*nnptfrx*(nnplan-1)
576  ELSE
577  WRITE(lu,*) 'ALMESH, UNEXPECTED ELEMENT FOR NELEB:',ielm
578  CALL plante(1)
579  stop
580  ENDIF
581 !
582  IF(reflevel.EQ.0) ALLOCATE(mesh%DIM1)
583  mesh%DIM1 = dimens(ielm0)
584  IF(reflevel.EQ.0) ALLOCATE(mesh%TYPELM)
585  mesh%TYPELM = ielm0
586  IF(reflevel.EQ.0) ALLOCATE(mesh%TYPELMBND)
587  mesh%TYPELMBND = typ_bnd_elem
588  IF(reflevel.EQ.0) ALLOCATE(mesh%NPOIN)
589  mesh%NPOIN = bief_nbpts(ielm1,mesh)
590  IF(reflevel.EQ.0) ALLOCATE(mesh%NPMAX)
591  mesh%NPMAX = bief_nbmpts(ielm1,mesh)
592  IF(reflevel.EQ.0) ALLOCATE(mesh%MXPTVS)
593  mesh%MXPTVS = mxptvs
594  IF(reflevel.EQ.0) ALLOCATE(mesh%MXELVS)
595  mesh%MXELVS = mxelvs
596 ! LV WILL BE RECOMPUTED LATER
597  IF(reflevel.EQ.0) ALLOCATE(mesh%LV)
598  mesh%LV = 1
599  IF(reflevel.EQ.0) ALLOCATE(mesh%NSEG)
600  mesh%NSEG = bief_nbseg(ielm1,mesh)
601 !
602 ! 2) ARRAYS OF INTEGERS
603 !
604 ! ALLOCATES IKLE AND KLEI (SAME SIZE, 2 INVERTED DIMENSIONS)
605 !
606  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%IKLE)
607  CALL bief_allvec(2,mesh%IKLE,'IKLE ',
608  & ielm0,bief_nbpel(ielm,mesh),1,mesh)
609  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%KLEI)
610  CALL bief_allvec(2,mesh%KLEI,'KLEI ',
611  & ielm0,bief_nbpel(ielm,mesh),1,mesh)
612 !
613 ! IFABOR
614 !
615  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%IFABOR)
616  CALL bief_allvec(2,mesh%IFABOR,'IFABOR',
617  & ielm0,bief_nbfel(ielm,mesh),1,mesh)
618 !
619 ! NELBOR
620 !
621 ! NELBOR & NULONE
622 ! IT IS NOW CHANGED TO VERTICAL BOUNDARY ELEMENT...
623 ! 2D NOT AFFECTED, 3D USAGE NELBO3(NPTFR,NETAGE) - NO. OF LAT. BD. ELEMENTS
624 !
625  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NELBOR)
626  CALL bief_allvec(2,mesh%NELBOR,'NELBOR',ielbor(ielm0,2),1,1,mesh)
627 !
628 ! NULONE
629 !
630 !
631 ! EXTRAORDINARILY STRANGE GEOMETRICALLY
632 ! IN 2D NUMBER OF BOUNDARY NODES IS EQUAL TO THE NUMBER OF BOUNDARY
633 ! ELEMENTS... IN 3D IT IS NOT THE CASE!
634 ! NULONE 3D IS USED INTERNALLY AS: NULONE(NPTFR,NETAGE,4)
635 ! "ASSOCIATES THE LOCAL BOUNDARY NUMBERING TO LOCAL 3D NUMBERING"
636 !
637  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NULONE)
638 !
639  CALL bief_allvec(2,mesh%NULONE,'NULONE',
640  & ielbor(ielm0,2),
641  & bief_nbpel(ielbor(ielm,2),mesh),1,mesh)
642 !
643 !!! NOTE : FOR THE TETRAHEDRONS, THIS IS NO LONGER THE CASE. WE READ THE
644 ! BOUNDARY CONNECTIVITY TABLE BEFORE INITIALISING THE NUMBER OF
645 ! BOUNDARY NODES AND ELEMENTS. THIS IS WELL DEFINED NOW.
646 !
647 ! IN 2D IT IS CALL BIEF_ALLVEC(2, MESH%NULONE, 'NULONE', 0, 2, 1,MESH)
648 ! WHICH, IN 2D ONLY IS EQUIVALENT TO
649 ! CALL BIEF_ALLVEC(2, MESH%NULONE, 'NULONE', 1, 2, 1,MESH)
650 ! IN 3D IT IS CALL BIEF_ALLVEC(2, MESH%NULONE, 'NULONE', 20, 4, 1,MESH)
651 !
652 !
653 ! KP1BOR
654 !
655  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%KP1BOR)
656  CALL bief_allvec(2,mesh%KP1BOR,'KP1BOR', ielbor(ielm,1),2,1,mesh)
657 !
658 ! NBOR: GLOBAL NUMBERS OF THE BOUNDARY NODES
659 ! ALLOCATES NBOR, IT WILL BE READ IN LECLIM
660 !
661  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NBOR)
662  CALL bief_allvec(2,mesh%NBOR,'NBOR ',ielbor(ielm,2),1,1,mesh)
663 !
664 ! IKLBOR: IKLE FOR THE SEGMENTS OR BOUNDARY SIDES
665 !
666  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%IKLBOR)
667 !
668  IF(ielm.EQ.41.OR.ielm.EQ.51) THEN
669 ! LATERAL BOUNDARY SIDES
670 ! SEE ININDS FOR ELEMENTS 20 AND 21
671 ! CALL BIEF_ALLVEC(2,MESH%IKLBOR,'IKLBOR',
672 ! & IELBOR(IELM0,2),
673 ! & BIEF_NBPEL(IELBOR(IELM1,2),MESH),1,MESH)
674  CALL bief_allvec(2,mesh%IKLBOR,'IKLBOR',
675  & mesh%NELEBX,
676  & bief_nbpel(ielbor(ielm1,2),mesh),0,mesh)
677  ELSEIF(ielm.EQ.11.OR.ielm.EQ.12.OR.ielm.EQ.13
678  & .OR.ielm.EQ.31) THEN
679  CALL bief_allvec(2,mesh%IKLBOR,'IKLBOR',
680  & ielbor(ielm0,1),
681  & bief_nbpel(ielbor(ielm ,2),mesh),1,mesh)
682  ELSE
683  WRITE(lu,*) 'ALMESH : UNKNOWN ELEMENT FOR IKLBOR:',ielm
684  CALL plante(1)
685  stop
686  ENDIF
687 !
688 ! IFANUM: NUMBER OF THE SIDE IN ADJACENT ELEMENT
689 !
690  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%IFANUM)
691  IF(cfg(1).EQ.2) THEN
692  CALL bief_allvec(2,mesh%IFANUM,'IFANUM',
693  & 3*bief_nbmpts(ielm0,mesh)
694  & + bief_nbmpts(01 ,mesh) ,1,0 ,mesh)
695  ELSEIF(cfg(1).NE.1.AND.cfg(1).NE.3) THEN
696  WRITE(lu,99) cfg(1)
697 99 FORMAT(1x,'ALMESH : UNKNOWN STORAGE:',1i6)
698  CALL plante(1)
699  stop
700  ELSE
701  CALL bief_allvec(2,mesh%IFANUM,'IFANUM', 0 , 1,0,mesh )
702  ENDIF
703 !
704 ! IKLEM1: INVERSE CONNECTIVITY TABLE FOR FRONTAL PRODUCT
705 ! LIMVOI: LIMITING NUMBER OF A GIVEN NUMBER OF NEIGHBOURS
706 !
707  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%IKLEM1)
708  CALL first_all_biefobj(mesh%LIMVOI)
709  IF(cfg(2).EQ.2) THEN
710 ! CALL BIEF_ALLVEC(2,MESH%IKLEM1,'IKLEM1',
711 ! NBMPTS(IELM1)*MXPTVS,4,0,MESH)
712 ! FOR OPTASS=3: SYM AND NOT SYM ARE DIFFERENT
713  CALL bief_allvec(2,mesh%IKLEM1,'IKLEM1',
714  & bief_nbmpts(ielm1,mesh)*mxptvs,8,0,mesh)
715 ! 11 IS HERE THE MAXIMUM MXPTVS PROGRAMMED IN OPASS
716  CALL bief_allvec(2,mesh%LIMVOI,'LIMVOI',11,2,0,mesh)
717  ELSEIF(cfg(2).NE.1) THEN
718  WRITE(lu,97) cfg(2)
719 97 FORMAT(1x,'ALMESH: UNKNOWN MATRIX-VECTOR PRODUCT:',1i6)
720  CALL plante(1)
721  stop
722  ELSE
723  CALL bief_allvec(2,mesh%IKLEM1,'IKLEM1',0,4,0,mesh)
724  CALL bief_allvec(2,mesh%LIMVOI,'LIMVOI',0,2,0,mesh)
725  ENDIF
726 !
727 ! INTEGER ARRAYS FOR SEGMENT-BASED STORAGE
728 !
729  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%GLOSEG)
730  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%ELTSEG)
731  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%ORISEG)
732  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%GLOSEGBOR)
733  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%ELTSEGBOR)
734  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%ORISEGBOR)
735 !
736  CALL bief_allvec(2,mesh%GLOSEG,'GLOSEG',
737  & bief_nbseg(ielm,mesh),2,0,mesh)
738  CALL bief_allvec(2,mesh%ELTSEG,'ELTSEG',
739  & bief_nbmpts(ielm0,mesh),
740  & bief_nbsegel(ielm,mesh),0,mesh)
741  CALL bief_allvec(2,mesh%ORISEG,'ORISEG',
742  & bief_nbmpts(ielm0,mesh),
743  & bief_nbsegel(ielm,mesh),0,mesh)
744  CALL bief_allvec(2,mesh%GLOSEGBOR,'GSGBOR',
745  & bief_nbseg(ielb1,mesh),2,0,mesh)
746  CALL bief_allvec(2,mesh%ELTSEGBOR,'ESGBOR',
747  & bief_nbmpts(ielb0,mesh),
748  & bief_nbsegel(ielb1,mesh),0,mesh)
749  CALL bief_allvec(2,mesh%ORISEGBOR,'OSGBOR',
750  & bief_nbmpts(ielb0,mesh),
751  & bief_nbsegel(ielb1,mesh),0,mesh)
752 !
753 ! INTEGER ARRAY FOR THE METHOD OF CHARACTERISTICS
754 ! THE STARTING ELEMENT (0 IF NOT IN THIS SUBDOMAIN)
755 !
756  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%ELTCAR)
757  CALL bief_allvec(2,mesh%ELTCAR,'ELTCAR',ielm,1,1,mesh)
758 !
759 ! INTEGER ARRAYS FOR PARALLEL MODE
760 !
761 ! KNOLG
762 ! NACHB
763 ! ISEG
764 ! INDPU
765 ! NHP
766 ! NHM
767 !
768  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%KNOLG)
769  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NACHB)
770  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%ISEG)
771  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%INDPU)
772  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NHP)
773  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NHM)
774  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%IFAPAR)
775  IF (reflevel.EQ.0) ALLOCATE(mesh%NB_NEIGHB)
776  IF (reflevel.EQ.0) ALLOCATE(mesh%NB_NEIGHB_SEG)
777 ! THERE ALLVEC IS IN PARINI
778  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NB_NEIGHB_PT)
779  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%LIST_SEND)
780  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NH_COM)
781  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NB_NEIGHB_PT_SEG)
782  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%LIST_SEND_SEG)
783  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NH_COM_SEG)
784 !
785 ! DATA STRUCTURE IN PARALLEL
786 !
787  IF(ncsize.GT.1) THEN
788 !
789  CALL bief_allvec(2,mesh%KNOLG,'KNOLG ',
790  & ielm1 ,1,1,mesh)
791  CALL bief_allvec(2,mesh%NACHB,'NACHB ',
792  & nbmaxnshare*nptir,1,0,mesh)
793  CALL bief_allvec(2,mesh%ISEG ,'ISEG ',
794  & ielbor(ielm1,1) ,1,1,mesh)
795  CALL bief_allvec(2,mesh%INDPU,'INDPU ',
796  & ielm1 ,1,1,mesh)
797  CALL bief_allvec(2,mesh%NHP ,'NHP ',
798  & nbmaxdshare*2*nptir,1,0,mesh)
799  CALL bief_allvec(2,mesh%NHM ,'NHM ',
800  & nbmaxdshare*2*nptir,1,0,mesh)
801  CALL bief_allvec(2,mesh%IFAPAR,'IFAPAR',10,6,1,mesh)
802 !
803  DO i=1,6*bief_nbpts(10,mesh)
804  mesh%IFAPAR%I(i)=0
805  ENDDO
806 !
807  ELSE
808 !
809  CALL bief_allvec(2,mesh%KNOLG ,'KNOLG ',0,1,0,mesh)
810  CALL bief_allvec(2,mesh%NACHB ,'NACHB ',0,1,0,mesh)
811  CALL bief_allvec(2,mesh%ISEG ,'ISEG ',0,1,0,mesh)
812  CALL bief_allvec(2,mesh%INDPU ,'INDPU ',0,1,0,mesh)
813  CALL bief_allvec(2,mesh%NHP ,'NHP ',0,1,0,mesh)
814  CALL bief_allvec(2,mesh%NHM ,'NHM ',0,1,0,mesh)
815  CALL bief_allvec(2,mesh%IFAPAR,'IFAPAR',0,1,0,mesh)
816 !
817  ENDIF
818 !
819 !-----------------------------------------------------------------------
820 !
821 ! FINITE VOLUMES
822 !
823  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%NUBO)
824  IF(equa(1:15).EQ.'SAINT-VENANT VF') THEN
825  CALL bief_allvec(2,mesh%NUBO,'NUBO ',2*nseg,1,0,mesh)
826  ELSE
827  CALL bief_allvec(2,mesh%NUBO,'NUBO ', 0,1,0,mesh)
828  ENDIF
829 !
830  IF (reflevel.EQ.0) CALL first_all_biefobj(mesh%JMI)
831  IF(equa(1:15).EQ.'SAINT-VENANT VF') THEN
832  CALL bief_allvec(2,mesh%JMI,'JMI ',nseg,1,0,mesh)
833  ELSE
834  CALL bief_allvec(2,mesh%JMI,'JMI ', 0,1,0,mesh)
835  ENDIF
836 !
837 !-----------------------------------------------------------------------
838 !
839 ! PART SKIPPED WHEN STARTING REFINEMENT
840 !
841  IF(reflevel.EQ.0) THEN
842 !
843 !-----------------------------------------------------------------------
844 !
845 ! FILLS ARRAYS IKLE, X AND Y (AND Z)
846 !
847 ! READS TO KNOLG ARRAY LOCAL TO GLOBAL NUMBERING IN PARALLEL
848  IF(ncsize.GT.1) THEN
849  CALL get_mesh_l2g_numbering(fformat,nfic,mesh%KNOLG%I,
850  & npoin,ierr)
851  CALL check_call(ierr,'ALMESH:GET_MESH_L2G_NUMBERING')
852  ENDIF
853 ! PRISMS CUT INTO TETRAHEDRONS
854 !
855  IF(ielm.EQ.51) THEN
856 !
857 ! NOTE : NO Z HERE, AS IELM.EQ.41, SEE NOTE BELOW
858  CALL read_mesh_coord(fformat,nfic,mesh%X%R,mesh%Y%R,npoin,
859  & projec,lati,longi)
860 !
861  CALL cpikle3(mesh%IKLE%I,ikles,nelem,nnelmax,npoin,nnplan,
862  & mesh%KNOLG%I)
863 !
864 ! PRISMS
865 !
866  ELSEIF(ielm.EQ.41) THEN
867 !
868 ! NOTE : WITH PRISMS Z IS COMPUTED WITH ZF AND H, OR
869 ! READ IN THE PREVIOUS COMPUTATION FILE, HENCE NO Z HERE
870  CALL read_mesh_coord(fformat,nfic,mesh%X%R,mesh%Y%R,npoin,
871  & projec,lati,longi)
872 !
873  CALL cpikle2(mesh%IKLE%I,mesh%KLEI%I,ikles,
874  & nelem,nnelmax,npoin,nnplan)
875 !
876 ! TRIANGLES OR TETRAHEDRONS
877 !
878  ELSEIF(ielm.EQ.11.OR.ielm.EQ.12.OR.ielm.EQ.13.OR.ielm.EQ.14
879  & .OR.ielm.EQ.31) THEN
880 !
881 ! IKLES(NDP,NELEM) COPIED INTO IKLE(NELMAX,NDP) AND KLEI(NDP,NELMAX)
882  DO i = 1,ndp
883  DO ielem = 1,nelem
884  mesh%IKLE%I((i-1)*nnelmax+ielem) = ikles((ielem-1)*ndp+i)
885  mesh%KLEI%I((ielem-1)*ndp+i) = ikles((ielem-1)*ndp+i)
886  ENDDO
887  ENDDO
888  IF(ielm.EQ.11.OR.ielm.EQ.12.OR.ielm.EQ.13.OR.ielm.EQ.14) THEN
889  CALL read_mesh_coord(fformat,nfic,mesh%X%R,mesh%Y%R,npoin,
890  & projec,lati,longi)
891  ELSEIF(ielm.EQ.31) THEN
892 ! TETRAHEDRONS: READS THE Z COORDINATE AFTER X AND Y
893  CALL read_mesh_coord(fformat,nfic,mesh%X%R,mesh%Y%R,npoin,
894  & projec,lati,longi,z=mesh%Z%R)
895  ENDIF
896 !
897  ELSE
898 !
899 ! OTHER ELEMENT TYPES
900 !
901  WRITE(lu,*) 'ALMESH: UNKNOWN ELEMENT:',ielm
902  CALL plante(1)
903  stop
904 !
905  ENDIF
906 !
907 ! NOW WE HAVE KNOLG IN 2D
908 ! FINDING THE NUMBER OF POINTS IN THE ORIGINAL MESH
909 ! = MAXIMUM OF ALL KNOLG OF ALL SUB-DOMAINS
910 !
911  npoin_glob=0
912  IF(ncsize.GT.1) THEN
913  DO i=1,npoin
914  npoin_glob=max(npoin_glob,mesh%KNOLG%I(i))
915  ENDDO
916  npoin_glob=p_max(npoin_glob)
917  ENDIF
918 !
919 ! COMPLEMENTS: ARRAYS X, Y FOR PRISMS AND TETRAHEDRONS
920 ! KNOLG FOR PRISMS AND TETRAHEDRONS
921 !
922  IF(ielm.EQ.41.OR.ielm.EQ.51) THEN
923  DO iplan = 2,nnplan
924  CALL ov_2( 'X=Y ' , mesh%X%R,iplan, mesh%X%R,1,
925  & mesh%X%R,1, 0.d0, nnpmax,npoin)
926  CALL ov_2( 'X=Y ' , mesh%Y%R,iplan, mesh%Y%R,1,
927  & mesh%Y%R,1, 0.d0, nnpmax,npoin)
928  ENDDO
929  IF(ncsize.GT.1) THEN
930  DO iplan = 2,nnplan
931  DO i=1,npoin
932  mesh%KNOLG%I(i+(iplan-1)*npoin)=
933  & mesh%KNOLG%I(i)+(iplan-1)*npoin_glob
934  ENDDO
935  ENDDO
936  ENDIF
937  ENDIF
938 !
939 !-----------------------------------------------------------------------
940 !
941 ! END OF PART SKIPPED WHEN STARTING REFINEMENT
942 !
943  ENDIF
944 !
945 !-----------------------------------------------------------------------
946 !
947 !
948 !-----------------------------------------------------------------------
949 !
950 ! DEALLOCATES TEMPORARY ARRAYS
951 !
952  DEALLOCATE(ikles)
953 !
954 !-----------------------------------------------------------------------
955 !
956  WRITE(lu,*) 'MESH: ',nom,' ALLOCATED'
957 !
958 !-----------------------------------------------------------------------
959 !
960  RETURN
961  END
subroutine mxptel(MXPTVS, MXELVS, IKLES, IELM, NPOIN, NELEM, NDP, IPOBO, LISTIN)
Definition: mxptel.f:7
integer function dimens(IELM)
Definition: dimens.f:7
integer function bief_nbpts(IELM, MESH)
Definition: bief_nbpts.f:7
integer function bief_nbpel(IELM, MESH)
Definition: bief_nbpel.f:7
integer function bief_nbfel(IELM, MESH)
Definition: bief_nbfel.f:7
subroutine bief_allvec(NAT, VEC, NOM, IELM, DIM2, STATUT, MESH)
Definition: bief_allvec.f:7
subroutine ov_2(OP, X, DIMX, Y, DIMY, Z, DIMZ, C, DIM1, NPOIN)
Definition: ov_2.f:7
subroutine bief_allmat(MAT, NOM, IELM1, IELM2, CFG, TYPDIA, TYPEXT, MESH)
Definition: bief_allmat.f:7
subroutine read_mesh_info(FFORMAT, NFIC, TITLE, NVAR, NPOIN, TYP_ELEM, NELEM, NPTFR, NPTIR, NDP, NPLAN, X_ORIG, Y_ORIG, TYP_BND_ELEM, NELEBD)
Definition: read_mesh_info.f:8
subroutine bief_ininds(NPOIN, NPTFR, NELEM, NPMAX, NPTFX, NELMAX, NPLAN, NSEGBOR, NDS, NELEB)
Definition: bief_ininds.f:7
subroutine segbor(NSEGBOR, IKLES, NELEM, NPOIN)
Definition: segbor.f:7
integer function ielbor(IELM, I)
Definition: ielbor.f:7
subroutine read_mesh_coord(FFORMAT, NFIC, X, Y, NPOIN, PROJECTION, LATI0, LONGI0, Z)
subroutine cpstvc(X, Y)
Definition: cpstvc.f:7
subroutine cpikle2(IKLE3, KLEI3, IKLES, NELEM2, NELMAX2, NPOIN2, NPLAN)
Definition: cpikle2.f:7
integer function bief_nbsegel(IELM, MESH)
Definition: bief_nbsegel.f:7
subroutine first_all_biefobj(OBJ)
subroutine almesh(MESH, NOM, IELM, SPHERI, CFG, FFORMAT, NFIC, EQUA, REFINE, NPLAN, NPMAX, NPTFRX, NELMAX, PROJECTION, LATI0, LONGI0, CONVERGENCE, RLEVEL)
Definition: almesh.f:8
subroutine read_mesh_conn(FFORMAT, NFIC, NPOIN, TYP_ELEM, NELEM, NDP, TYP_BND_ELEM, NELEBD, IKLE, IPOBO)
Definition: read_mesh_conn.f:8
integer function bief_nbmpts(IELM, MESH)
Definition: bief_nbmpts.f:7
subroutine cpikle3(IKLE3, IKLES, NELEM2, NELMAX2, NPOIN2, NPLAN, KNOLG)
Definition: cpikle3.f:7
subroutine get_mesh_l2g_numbering(FFORMAT, FID, KNOLG, NPOIN, IERR)
integer function bief_nbseg(IELM, MESH)
Definition: bief_nbseg.f:7
Definition: bief.f:3