Skip to content

GitLab

Commit be532f60 authored by mozul's avatar mozul
Browse files

Merge commit '0f711c6a' 

Merge dev sources
parents 95c16cb1 0f711c6a
Expand all Hide whitespace changes
Inline Side-by-side
src/ChiPy/src/macro/chipy.py
View file @ be532f60
......@@ -1018,7 +1018,7 @@ def WriteLastVlocRloc():
SPPLx_WriteLastVlocRloc()
SPSPx_WriteLastVlocRloc()
def WriteOutDof(nsteps):
def WriteOutDof(nsteps=1):
"""Write degrees of freedom of current time step in ascii file (OUTBOX/DOF.x.OUT where 'x' is the rank of the file)
if current time step number is multiple of nsteps.
"""
......@@ -1051,13 +1051,13 @@ def OpenDisplayFiles(restart=1,ref_reac=0.):
if is_vtk_display :
wd = overall_GetWorkingDirectory()
fim = startCollection(wdf,os.path.join(wd,'DISPLAY','mecafe.pvd'))
fin = startCollection(wdf,os.path.join(wd,'DISPLAY','therfe.pvd'))
fip = startCollection(wdf,os.path.join(wd,'DISPLAY','porofe.pvd'))
fiu = startCollection(wdf,os.path.join(wd,'DISPLAY','multife.pvd'))
fit = startCollection(wdf,os.path.join(wd,'DISPLAY','tacts.pvd'))
fii = startCollection(wdf,os.path.join(wd,'DISPLAY','inters.pvd'))
fir = startCollection(wdf,os.path.join(wd,'DISPLAY','rigids.pvd'))
fim = startCollection(os.path.join(wd,'DISPLAY','mecafe.pvd'),wdf)
fin = startCollection(os.path.join(wd,'DISPLAY','therfe.pvd'),wdf)
fip = startCollection(os.path.join(wd,'DISPLAY','porofe.pvd'),wdf)
fiu = startCollection(os.path.join(wd,'DISPLAY','multife.pvd'),wdf)
fit = startCollection(os.path.join(wd,'DISPLAY','tacts.pvd'),wdf)
fii = startCollection(os.path.join(wd,'DISPLAY','inters.pvd'),wdf)
fir = startCollection(os.path.join(wd,'DISPLAY','rigids.pvd'),wdf)
if DIMENSION == 2:
tact_names = ['DISKx', 'DISKx', 'DISPx', 'JONCx', 'POLYG', 'PT2Dx', 'xKSID', 'xPSID']
inter_names = ['CLALp', 'CLJCx', 'DKALp', 'DKDKL', 'DKDKx', 'DKDPx', 'DKJCx', 'DKKDx',
......@@ -1179,7 +1179,7 @@ def ClosePostproFiles():
elif DIMENSION == 3:
postpro_3D_ClosePostproFiles()
def WriteOutGPV(nsteps):
def WriteOutGPV(nsteps=1):
"""Write current Gauss points values of all meshes (OUTBOX/GPV.x.OUT)
if current time step number is a multiple of nsteps.
"""
......@@ -1196,7 +1196,7 @@ def WriteOutMpValues():
elif DIMENSION == 3:
mp_solver_3D_WriteOutMpValues()
def WriteOutRnod(nsteps):
def WriteOutRnod(nsteps=1):
"""Write current values of reaction for each bodies if current time step number
is a multiple of nsteps.
"""
......@@ -1210,7 +1210,7 @@ def WriteOutRnod(nsteps):
therMAILx_WriteOutRnod()
poroMAILx_WriteOutRnod()
def WriteOutVlocRloc(nsteps):
def WriteOutVlocRloc(nsteps=1):
"""Write values of every interactions if current time step number is a multiple of nsteps.
"""
global DIMENSION
......
src/ChiPy/src/macro/ddm_utils.py
View file @ be532f60
......@@ -4,7 +4,11 @@ import os
import subprocess
import shutil
import glob
import vtk
try:
import vtk
is_vtk_display = True
except ImportError:
is_vtk_display = False
import lmgc90
......@@ -79,6 +83,9 @@ def concatenate_DISPLAY():
"""
import glob
if not is_vtk_dispaly:
return
working_dir = lmgc90.overall_GetWorkingDirectory()
# Calcul du nombre de pas de temps pour lesquels des sorties
......@@ -169,5 +176,4 @@ def concatenate_DISPLAY():
my_output.write("</Collection>")
my_output.write("</VTKFile>")
my_output.close()
\ No newline at end of file
src/ChiPy/src/macro/vtk_display.py
View file @ be532f60
......@@ -530,7 +530,7 @@ def InitInter(dict,name):
dict[name]=vals
def pushInters2D(inters_dict,Append,lr,normal_orient,kw):
def pushInters2D(inters_dict,Append,xxAppend,lr,normal_orient,kw):
nb_write=0
......@@ -560,6 +560,8 @@ def pushInters2D(inters_dict,Append,lr,normal_orient,kw):
# pour eviter de diviser par 0. apres
if seuil_fn == 0.: seuil_fn=1.
if ref_reac != 0.: seuil_fn = ref_reac
if len(kw) != 0:
for name,dico in kw.items():
#print '-> ',name
......@@ -574,9 +576,12 @@ def pushInters2D(inters_dict,Append,lr,normal_orient,kw):
# on collecte
points = vtkPoints()
cells = vtkCellArray()
Ids = vtkFloatArray()
Ids.SetNumberOfComponents(1)
Ids.SetName("Ids")
Ids.InsertNextTuple1(float(iinter))
Rn = vtkFloatArray()
Rn.SetNumberOfComponents(1)
Rn.SetName("Rn")
......@@ -603,8 +608,6 @@ def pushInters2D(inters_dict,Append,lr,normal_orient,kw):
Gap.SetName("Gap")
Gap.InsertNextTuple1(g);
if ref_reac != 0.: seuil_fn = ref_reac
if normal_orient:
vertex = buildFrite(coorx,coory,nx,ny,fn/seuil_fn,lr)
else:
......@@ -621,7 +624,7 @@ def pushInters2D(inters_dict,Append,lr,normal_orient,kw):
num_points = 0
cell_ids=[]
for x,y in vertex:
points.InsertPoint(num_points,x,y,0)
points.InsertPoint(num_points,x,y,0.)
cell_ids.append(num_points)
num_points+=1
......@@ -632,7 +635,7 @@ def pushInters2D(inters_dict,Append,lr,normal_orient,kw):
cell.GetPointIds().SetId(i,j)
cells.InsertNextCell(cell)
Ids.InsertNextTuple1(float(iinter))
# ca pousse dans un polydata
polydata = vtkPolyData()
......@@ -659,9 +662,134 @@ def pushInters2D(inters_dict,Append,lr,normal_orient,kw):
add_data_to_append(Append, polydata)
###
ptc = vtkPoints()
ptc.InsertNextPoint(coorx,coory,0.)
Id = vtkFloatArray()
Id.SetNumberOfComponents(1)
Id.SetName("Ids")
Id.InsertNextTuple1(float(iinter))
N = vtkFloatArray()
N.SetNumberOfComponents(3)
N.SetName("N")
N.InsertNextTuple3(nx/lr,ny/lr, 0.)
R = vtkFloatArray()
R.SetNumberOfComponents(3)
R.SetName("R")
R.InsertNextTuple3((fn*nx + ft*ny)/(seuil_fn*lr), (fn*ny - ft*nx)/(seuil_fn*lr), 0.)
cell = vtkVertex()
cell.GetPointIds().SetId(0,0)
ug= vtkUnstructuredGrid()
ug.SetPoints(ptc)
ug.InsertNextCell(cell.GetCellType(),cell.GetPointIds())
ug.GetCellData().AddArray(Id)
ug.GetPointData().AddArray(N)
ug.GetPointData().AddArray(R)
add_data_to_append(xxAppend, ug)
if nb_write == 0:
Ids = vtkFloatArray()
Ids.SetNumberOfComponents(1)
Ids.SetName("Ids")
Rn = vtkFloatArray()
Rn.SetNumberOfComponents(1)
Rn.SetName("Rn")
Rt = vtkFloatArray()
Rt.SetNumberOfComponents(1)
Rt.SetName("Rt")
Strong = vtkUnsignedIntArray()
Strong.SetNumberOfComponents(1)
Strong.SetName("Strong")
Gap = vtkFloatArray()
Gap.SetNumberOfComponents(1)
Gap.SetName("Gap")
vertex = buildFrite(0.,0.,1.,0.,0.,1e-14)
# on pose les 4 points de la frite
cells = vtkCellArray()
points = vtkPoints()
num_points = 0
cell_ids=[]
for x,y in vertex:
points.InsertPoint(num_points,x,y,0.)
cell_ids.append(num_points)
num_points+=1
cell = vtkPolygon()
cell.GetPointIds().SetNumberOfIds(len(cell_ids))
for i,j in enumerate(cell_ids):
cell.GetPointIds().SetId(i,j)
cells.InsertNextCell(cell)
Ids.InsertNextTuple1(float(0))
Rn.InsertNextTuple1(0.);
Rt.InsertNextTuple1(0.);
Strong.InsertNextTuple1(0);
Gap.InsertNextTuple1(0.);
# ca pousse dans un polydata
polydata = vtkPolyData()
polydata.SetPoints(points)
polydata.SetPolys(cells)
polydata.GetCellData().AddArray(Ids)
polydata.GetCellData().AddArray(Rn)
polydata.GetCellData().AddArray(Rt)
polydata.GetCellData().AddArray(Strong)
polydata.GetCellData().AddArray(Gap)
add_data_to_append(Append, polydata)
###
ptc = vtkPoints()
ptc.InsertNextPoint(0.,0.,0.)
Id = vtkFloatArray()
Id.SetNumberOfComponents(1)
Id.SetName("Ids")
Id.InsertNextTuple1(float(0.))
N = vtkFloatArray()
N.SetNumberOfComponents(3)
N.SetName("N")
N.InsertNextTuple3(0.,0.,0.)
R = vtkFloatArray()
R.SetNumberOfComponents(3)
R.SetName("R")
R.InsertNextTuple3(0.,0.,0.)
cell = vtkVertex()
cell.GetPointIds().SetId(0,0)
ug= vtkUnstructuredGrid()
ug.SetPoints(ptc)
ug.InsertNextCell(cell.GetCellType(),cell.GetPointIds())
ug.GetCellData().AddArray(Id)
ug.GetPointData().AddArray(N)
ug.GetPointData().AddArray(R)
add_data_to_append(xxAppend, ug)
nb_write=1
return nb_write
def pushInters3D(inters_dict,Append,lr,normal_orient,kw):
def pushInters3D(inters_dict,Append,xxAppend,lr,normal_orient,kw):
nb_write=0
......@@ -695,6 +823,8 @@ def pushInters3D(inters_dict,Append,lr,normal_orient,kw):
# pour eviter de diviser par 0. apres
if seuil_fn == 0.: seuil_fn=1.
if ref_reac != 0.: seuil_fn = ref_reac
if len(kw) != 0:
for name,dico in kw.items():
#print '-> ',name
......@@ -731,8 +861,6 @@ def pushInters3D(inters_dict,Append,lr,normal_orient,kw):
Gap.SetNumberOfComponents(1)
Gap.SetName("Gap")
if ref_reac != 0.: seuil_fn = ref_reac
if normal_orient:
vertex = buildFrite3D(coorx,coory,coorz,nx,ny,nz,tx,ty,tz,fn/seuil_fn,lr)
else:
......@@ -801,6 +929,131 @@ def pushInters3D(inters_dict,Append,lr,normal_orient,kw):
add_data_to_append(Append, polydata)
###
ptc = vtkPoints()
ptc.InsertNextPoint(coorx,coory,coorz)
Id = vtkFloatArray()
Id.SetNumberOfComponents(1)
Id.SetName("Ids")
Id.InsertNextTuple1(float(iinter))
N = vtkFloatArray()
N.SetNumberOfComponents(3)
N.SetName("N")
N.InsertNextTuple3(nx/lr,ny/lr,nz/lr)
#N.InsertNextTuple3(nx,ny,nz)
R = vtkFloatArray()
R.SetNumberOfComponents(3)
R.SetName("R")
R.InsertNextTuple3((fn*nx + ft*tx)/(seuil_fn*lr), (fn*ny + ft*ty)/(seuil_fn*lr), (fn*nz + ft*tz)/(seuil_fn*lr))
#R.InsertNextTuple3((fn*nx + ft*tx), (fn*ny + ft*ty), (fn*nz + ft*tz))
cell = vtkVertex()
cell.GetPointIds().SetId(0,0)
ug= vtkUnstructuredGrid()
ug.SetPoints(ptc)
ug.InsertNextCell(cell.GetCellType(),cell.GetPointIds())
ug.GetCellData().AddArray(Id)
ug.GetPointData().AddArray(N)
ug.GetPointData().AddArray(R)
add_data_to_append(xxAppend, ug)
if nb_write == 0:
Ids = vtkFloatArray()
Ids.SetNumberOfComponents(1)
Ids.SetName("Ids")
Rn = vtkFloatArray()
Rn.SetNumberOfComponents(1)
Rn.SetName("Rn")
Rt = vtkFloatArray()
Rt.SetNumberOfComponents(1)
Rt.SetName("Rt")
Strong = vtkUnsignedIntArray()
Strong.SetNumberOfComponents(1)
Strong.SetName("Strong")
Gap = vtkFloatArray()
Gap.SetNumberOfComponents(1)
Gap.SetName("Gap")
vertex = buildFrite3D(0.,0.,0.,1.,0.,0.,0.,1.,0.,0.,1e-14)
# on pose les 8 points de la frite
num_points = 0
points = vtkPoints()
for x,y,z in vertex:
points.InsertPoint(num_points,x,y,z)
num_points+=1
# on va creer un patch avec les 6 faces de la frite (un hexa ca ne marche pas)
# on doit passer a chaque frite un tableau de valeur par face doit contenir les valeurs
cells = vtkCellArray()
for j in range(6):
cell = vtkQuad()
for k in range(4):
cell.GetPointIds().SetId(k,frite_connec[j,k])
cells.InsertNextCell(cell)
Ids.InsertNextTuple1(float(0))
Rn.InsertNextTuple1(0.);
Rt.InsertNextTuple1(0.);
Strong.InsertNextTuple1(0);
Gap.InsertNextTuple1(0.);
# ca pousse dans un polydata
polydata = vtkPolyData()
polydata.SetPoints(points)
polydata.SetPolys(cells)
polydata.GetCellData().AddArray(Ids)
polydata.GetCellData().AddArray(Rn)
polydata.GetCellData().AddArray(Rt)
polydata.GetCellData().AddArray(Strong)
polydata.GetCellData().AddArray(Gap)
add_data_to_append(Append, polydata)
###
ptc = vtkPoints()
ptc.InsertNextPoint(0.,0.,0.)
Id = vtkFloatArray()
Id.SetNumberOfComponents(1)
Id.SetName("Ids")
Id.InsertNextTuple1(float(0.))
N = vtkFloatArray()
N.SetNumberOfComponents(3)
N.SetName("N")
N.InsertNextTuple3(0.,0.,0.)
R = vtkFloatArray()
R.SetNumberOfComponents(3)
R.SetName("R")
R.InsertNextTuple3(0.,0.,0.)
cell = vtkVertex()
cell.GetPointIds().SetId(0,0)
ug= vtkUnstructuredGrid()
ug.SetPoints(ptc)
ug.InsertNextCell(cell.GetCellType(),cell.GetPointIds())
ug.GetCellData().AddArray(Id)
ug.GetPointData().AddArray(N)
ug.GetPointData().AddArray(R)
add_data_to_append(xxAppend, ug)
nb_write=1
return nb_write
def pushNewInters2D(Append,lr,kw):
......@@ -1049,12 +1302,20 @@ def writeIntersToVTK(fichier,fid,inters_dict,lr,normal_orient=True,**kw):
AppendAll=vtkAppendPolyData()
# < xxx
vtuFile = vtkXMLUnstructuredGridWriter()
xxfichier=fichier.replace('inters','ptc')
xxfichier=xxfichier.replace('vtp','vtu')
vtuFile.SetFileName(xxfichier)
xxAppendAll=vtkAppendFilter()
# xxx >
if inters_dict.keys()[0] in ['CDCDx', 'CDPLx', 'CSASp', 'CSPRx', 'PRPLx', 'PRPRx', 'PTPT3', 'SPCDx', 'SPDCx', 'SPPLx', 'SPSPx']:
nb = pushInters3D(inters_dict,AppendAll,lr,normal_orient,kw)
nb = pushInters3D(inters_dict,AppendAll,xxAppendAll,lr,normal_orient,kw)
else:
nb = pushInters2D(inters_dict,AppendAll,lr,normal_orient,kw)
nb = pushInters2D(inters_dict,AppendAll,xxAppendAll,lr,normal_orient,kw)
# si on a bien ecrit qq ch
......@@ -1066,6 +1327,11 @@ def writeIntersToVTK(fichier,fid,inters_dict,lr,normal_orient=True,**kw):
impr='<DataSet timestep="%s" group="" part="0" file="%s"/>\n' % (time,'./'+os.path.basename(fichier))
fid.write(impr)
# xxx
set_input_of_file(vtuFile, xxAppendAll)
vtuFile.Write()
def writeNewIntersToVTK(f_path,f_name,fid,lr,dim,**kw):
vtpFile = vtkXMLPolyDataWriter()
......@@ -2298,7 +2564,7 @@ def writeMultifeToVTK(fichier,fid,dim,**kw):
#######################################
# ecriture collections
def startCollection(restart,fichier):
def startCollection(fichier,restart=0):
if restart > 0 :
# taken from stackoverflow :
# smart way to remove last ltd lines from a file
......
src/Core/src/contact_2D/mod_CLALp.f90
View file @ be532f60
......@@ -433,48 +433,6 @@ contains
end select
end subroutine write_xxx_Vloc_Rloc_CLALp
!!$!!!--------------------------------------------------------------
!!$ SUBROUTINE display_contact_energy_CLALp
!!$
!!$ IMPLICIT NONE
!!$
!!$ REAL(kind=8) :: NRJ,NRJ_t,NRJ_n
!!$
!!$!!!mr? Mais ou est le close????
!!$
!!$ IF (.NOT. is_open_energy) THEN
!!$ OPEN(56)
!!$ is_open_energy=.TRUE.
!!$ NRJ=0.d0
!!$ ENDIF
!!$
!!$ CALL compute_energy_increment(NRJ_t,NRJ_n)
!!$
!!$ NRJ = NRJ + (NRJ_t + NRJ_n)
!!$
!!$ WRITE (56,*) TPS,NRJ
!!$
!!$ END SUBROUTINE display_contact_energy_CLALp
!!$!!!--------------------------------------------------------------
!!$ SUBROUTINE display_contact_internal_CLALp
!!$
!!$ IMPLICIT NONE
!!$
!!$ INTEGER :: ik
!!$!!!mr? Tiens le copain du dessus!!!
!!$ ik = 1
!!$
!!$ IF (.NOT. is_open_internal) THEN
!!$ OPEN(55)
!!$ is_open_internal=.TRUE.
!!$ END IF
!!$
!!$ WRITE(55,'(6(1x,D14.6))') TPS,&
!!$ this(ik)%internal(2),this(ik)%internal(3),this(ik)%internal(4), &
!!$ this(ik)%rlt/H,this(ik)%rln/H
!!$
!!$ END SUBROUTINE display_contact_internal_CLALp
!!$!!!------------------------------------------------------------------------
!!!------------------------------------------------------------------------
subroutine set_nonsymmetric_detection_CLALp()
......@@ -1807,7 +1765,8 @@ contains
allocate(violation(nb_CLALp),stat=errare)
end subroutine compute_contact_CLALp
!------------------------------------------------------------------------
!------------------------------------------------------------------------
!> \brief Compute one contact from a CLALp rough
subroutine compute_one_contact_CLALp(i4_input, i4_output, r8_output, to_keep)
implicit none
......@@ -2043,7 +2002,6 @@ contains
end if
end subroutine compute_one_contact_CLALp
!------------------------------------------------------------------------
!------------------------------------------------------------------------
subroutine display_prox_tactors_CLALp
......@@ -3490,7 +3448,9 @@ end subroutine print_info_CLALp
end function get_all_internal_CLALp
!------------------------------------------------------------------------
!------------------------------------------------------------------------
!! ! rm : functions to test itrHdl
!! !> \brief Get a rough interaction
!! subroutine get_rough_CLALp(icdan, icdtac, iantac, isee, periodic, i4)
......
src/Core/src/contact_2D/mod_CLJCx.f90
View file @ be532f60
......@@ -1944,7 +1944,7 @@ END FUNCTION CHECK_CLJCx
get_iantac_CLJCx = this(icdan)%iantac
end function
!------------------------------------------------------------------------
subroutine clean_memory_CLJCx
implicit none
integer(kind=4) :: i
......
src/Core/src/contact_3D/mod_PRPLx.f90
View file @ be532f60
This diff is collapsed.
src/Core/src/contact_3D/mod_PRPRx.f90
View file @ be532f60
......@@ -296,7 +296,7 @@ MODULE PRPRx
! rank in this array for adjacent contact iadj
! to candidate contactor icdtac.
 
INTEGER :: cdtype ! model type =1 RBDY3, =2 MAILx
INTEGER :: cdtype ! model type =i_rbdy3 RBDY3, =i_mecaM MAILx , =i_mbs3 MBS3D
 
INTEGER :: cdbdy ! serial number of candidate body
 
......@@ -304,7 +304,7 @@ MODULE PRPRx
! By definition verlt(icdtac)%cdtac =icdtac;
 
INTEGER,DIMENSION(:) ,POINTER :: antype ! verlt(icdtac)%antype(iadj):
! model type =1 RBDY3, =2 MAILx of antagoniste
! model type = i_rbdy3 RBDY3, =i_mecaM MAILx, =i_mbs3 MBS3D of antagoniste
INTEGER,DIMENSION(:) ,POINTER :: anbdy ! verlt(icdtac)%anbdy(iadj):
! serial number of antagonist body for adjacent contactor iadj.
 
......@@ -836,8 +836,10 @@ SUBROUTINE binary_write_xxx_Vloc_Rloc_PRPRx(io_unit)
 
cdbdy = 'RBDY3'
if (polyr2bdyty(3,icdtac) == i_mecaM) cdbdy = 'MAILx'
if (polyr2bdyty(3,icdtac) == i_mbs3) cdbdy = 'MBS3D'
anbdy = 'RBDY3'
if (polyr2bdyty(3,iantac) == i_mecaM) anbdy = 'MAILx'
if (polyr2bdyty(3,iantac) == i_mbs3) anbdy = 'MBS3D'
!mr must be defined during the detection
icdver = 0
......@@ -902,7 +904,9 @@ SUBROUTINE display_prox_tactors_PRPRx
 
nb_POLYR=get_nb_POLYR()
 
DO icdtac=1,nb_POLYR
DO icdtac=1,nb_POLYR
DO iadj=1,nb_adj(icdtac)
icdan = adjac(icdtac)%icdan(iadj)
icdbdy = this(icdan)%icdbdy
......@@ -912,8 +916,11 @@ SUBROUTINE display_prox_tactors_PRPRx
 
cdtype='RBDY3'
if (polyr2bdyty(3,this(icdan)%icdtac) == i_mecaM) cdtype='MAILx'
if (polyr2bdyty(3,this(icdan)%icdtac) == i_mbs3) cdtype='MBS3D'
antype='RBDY3'
if (polyr2bdyty(3,this(icdan)%iantac) == i_mecaM) antype='MAILx'
if (polyr2bdyty(3,this(icdan)%iantac) == i_mbs3) antype='MBS3D'
 
WRITE(*,'(A1)')' '
WRITE(*,'(A6,2X,I5)')'$icdan',icdan
......@@ -1123,9 +1130,11 @@ SUBROUTINE stock_rloc_PRPRx
iadj = this(icdan)%iadj ! serial adjacent number of pair contactor
! adjacent to candidate contactor for contact icdan
verlt(icdtac)%icdan(iadj) = icdan