TABLE OF CONTENTS
UTI/PRE/init_mesh [ Modules ]
NOM
init_mesh (msh)
DESCRIPTION
Calcul les donnees geometriques (surfaces, volumes, normales)
a partir de donnees elementaires (liste des points et connectivite faces-cellules)
Calculate geometric data (surfaces, volumes, normals)
from basic data (list of points and face-cell connectivity)
ENTREES / INPUT
msh: le maillage / the mesh
SORTIES / OUTPUT
msh: le maillage / the mesh
SOURCE
! Ce logiciel est regi par la licence [CeCILL-B] ! This software is governed by the [CeCILL-B] license !=========================== DEBUT DES DECLARATIONS ==================== !.1----- Implicit, Use use geo_typ use geo, only : triangle IMPLICIT NONE !.2----- Declaration type(mesh), intent(inout) :: msh !! Maillage / Mesh real(kind=kind(0.d0)) :: x1,y1,z1,x2,y2,z2,x3,y3,z3,x4,y4,z4 real(kind=kind(0.d0)) :: xn1,yn1,zn1,surf1,xm1,ym1,zm1,xn2,yn2,zn2,surf2,xm2,ym2,zm2 real(kind=kind(0.d0)) :: xn,yn,zn,surf,xm,ym,zm,d13,d24 integer :: i1,i2,i3,i4,nf,iv,iv1,iv2 !=========================== DEBUT DU CODE EXECUTABLE ================== DO nf=1,msh%nb_face ! on decoupe chaque face en 2 triangles ! coplanaires !!!!! i1=msh%list_face(nf)%vertex(1) ! n° du sommet 1 i2=msh%list_face(nf)%vertex(2) i3=msh%list_face(nf)%vertex(3) x1=msh%list_vertex(i1)%x y1=msh%list_vertex(i1)%y z1=msh%list_vertex(i1)%z x2=msh%list_vertex(i2)%x y2=msh%list_vertex(i2)%y z2=msh%list_vertex(i2)%z x3=msh%list_vertex(i3)%x y3=msh%list_vertex(i3)%y z3=msh%list_vertex(i3)%z IF (msh%list_face(nf)%nbvertex.eq.4) THEN i4=msh%list_face(nf)%vertex(4) x4=msh%list_vertex(i4)%x y4=msh%list_vertex(i4)%y z4=msh%list_vertex(i4)%z d13=(x3-x1)**2+(y3-y1)**2+(z3-z1)**2 d24=(x4-x2)**2+(y4-y2)**2+(z4-z2)**2 IF (d13.le.d24) THEN CALL triangle(x1,y1,z1,x2,y2,z2,x3,y3,z3,xn1,yn1,zn1,surf1,xm1,ym1,zm1) CALL triangle(x1,y1,z1,x3,y3,z3,x4,y4,z4,xn2,yn2,zn2,surf2,xm2,ym2,zm2) ELSE CALL triangle(x1,y1,z1,x2,y2,z2,x4,y4,z4,xn1,yn1,zn1,surf1,xm1,ym1,zm1) CALL triangle(x2,y2,z2,x3,y3,z3,x4,y4,z4,xn2,yn2,zn2,surf2,xm2,ym2,zm2) END IF surf=surf1+surf2 IF (surf.ne.0) THEN xn=(surf1*xn1+surf2*xn2)/surf yn=(surf1*yn1+surf2*yn2)/surf zn=(surf1*zn1+surf2*zn2)/surf xm=(surf1*xm1+surf2*xm2)/surf ym=(surf1*ym1+surf2*ym2)/surf zm=(surf1*zm1+surf2*zm2)/surf ELSE xn=-1 yn=0 zn=0 xm=0 ym=0 zm=0 END IF ELSE CALL triangle(x1,y1,z1,x2,y2,z2,x3,y3,z3,xn,yn,zn,surf,xm,ym,zm) END IF msh%list_face(nf)%center%x=xm msh%list_face(nf)%center%y=ym msh%list_face(nf)%center%z=zm msh%list_face(nf)%surf=surf msh%list_face(nf)%vnorm%v(1)=xn msh%list_face(nf)%vnorm%v(2)=yn msh%list_face(nf)%vnorm%v(3)=zn END DO ! Calcul des rapports Volume/surface DO iv=1,msh%nb_cell IF(msh%list_cell(iv)%idom.eq.msh%numdom)msh%list_cell(iv)%h=0.d0 end DO DO nf=1,msh%nb_face iv1=msh%list_face(nf)%ic1 iv2=msh%list_face(nf)%ic2 IF (iv1.gt.0.and.iv2.ne.-999) THEN !!!!!!!!!!!!! ON ESSAYE !!!!!!!!!!! IF(msh%list_cell(iv1)%idom.eq.msh%numdom)msh%list_cell(iv1)%h=msh%list_cell(iv1)%h+msh%list_face(nf)%surf END IF IF (iv2.gt.0) THEN IF(msh%list_cell(iv2)%idom.eq.msh%numdom)msh%list_cell(iv2)%h=msh%list_cell(iv2)%h+msh%list_face(nf)%surf END IF END DO msh%small=1.d20 msh%big =0.d0 DO iv=1,msh%nb_cell IF(msh%list_cell(iv)%idom.eq.msh%numdom)then msh%list_cell(iv)%h=msh%list_cell(iv)%vol/msh%list_cell(iv)%h if(msh%list_cell(iv)%h.gt.msh%big) msh%big =msh%list_cell(iv)%h if(msh%list_cell(iv)%h.lt.msh%small)msh%small=msh%list_cell(iv)%h endif END DO write( *,1000)msh%numdom,msh%small write( *,1010)msh%numdom,msh%small 1000 format('INIT_MESH-INFO Domain ',i3,' size of smallest cell :',e13.5) 1010 format('INIT_MESH-INFO Domain ',i3,' size of largest cell :',e13.5) msh%list_face(:)%flevel=0 msh%list_cell(:)%clevel=0 msh%nb_level_cell=int(dlog(msh%big*1.05/msh%small)/dlog(2.0d0)) msh%nb_level_face=msh%nb_level_cell !=========================== FIN DE LA ROUTINE ==================== END SUBROUTINE INIT_MESH