← Back to team overview

yade-users team mailing list archive

[Question #708272]: i.geom.penetrationDepth = nan value when saving ViscElCap interactions

 

New question #708272 on Yade:
https://answers.launchpad.net/yade/+question/708272

Hello, I need to save then reload simulation with ViscElCap interaction inside but when I reload it the penetration Depth of ViscElCap interaction is equal to nan. This  nan value contaminates position and velocity values making disappear particles.

I softly fix that by erasing and recreate each ViscElCap interaction, but obviously it didn't conserve friction history which is important to me.

Maybe someone know why and how to fix this ?

I'm using yade 2022.01a 

You can test that with : 

#!/usr/bin/env python3
# -*- coding: utf-8 -*-


#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# import
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

import sys
import numpy as np
from yade import plot, qt


#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#DEF
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

s_name = sys.argv[0] [0:-3]


#sim param --------------------------------------------------------------------

p_speed = 0.01

tc_frac = 0.0002

#Mechanics --------------------------------------------------------------------

#Value in order :

#mean radius

R1 = 19 * 1e-6
R2 = 27.5 * 1e-6 #35 #32.5 #27.5

#radius = 2 * R1 * R2 /(R1 + R2)
radius = 0.0002 #0.002381 #2e-4

#Law {NoCoh,LiCoh,SoCoh}
Law = 'ViscElCap' #CohFrictMat

#particles
p_density = 2650 #2650
f_density = 1000

#Modulus
E_p = 1e8 #1e8

#poisson
P_p = 0.2

#friction angle
fa_p =  0.5

#restit
rcn_p = 0.7
rct_p = 0.7

#cohesion ---------------------------------------------------------------------

#contact time
tc = 0.001

#surface tension
Gamma = 72.8 * 1e-3 #27 #24 # 28

#Wetting angle [°]
Theta = 10  #0

#Liquid bridge volume def

#[%] of the Void volume filled with water
water_content = 0.1 #from [0.0001 ; 0.1]

Nl = 5.00481359595265 #6 #to precise 

#VB = 36 * 1e-19 #2 #12 #36
#VB = 2 * water_content * Vviod / Nl
VB = water_content * np.pi * p_density * radius**3 / ( 3 * Nl * f_density )


#seed type=int (no seed = -1)
seed_gen = 1

#Liquid Coh
CapModel = ['Willett_numeric',
            'Willett_analytic',
            'Weigert',
            'Rabinovich',
            'Lambert',
            'Soulie',]


SolidModel = ['']


#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#Pre compute
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

#create result folder if no existing

#print
print()

print('Radius : ' + str(radius) )
print('Bridge volume : ' + str(VB) )

print('\n\n\n\n\n\n\n\n\n\n')

#r.seed(1)


#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#Material def
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

#------------------------------------------------------------------------------
#Cohesive (Liquid)
#------------------------------------------------------------------------------

for i in range(len(CapModel)) :
    O.materials.append(ViscElCapMat(label = CapModel[i], frictionAngle = fa_p, mR = 0.1, mRtype = 1, en = rcn_p, et = rcn_p, tc = tc, Capillar = True, gamma = Gamma, Vb = VB, theta = 10, CapillarType = CapModel[i] ))

    
nbCm = len(CapModel)
Models = CapModel


#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#Engine def
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

#Engine def -------------------------------------------------------------------

O.engines = [ForceResetter(),
             InsertionSortCollider([Bo1_Sphere_Aabb(), Bo1_Facet_Aabb()]),
             
             #Collision detection
             InteractionLoop([Ig2_Sphere_Sphere_ScGeom(), Ig2_Facet_Sphere_ScGeom()],
                             
                             
                             [Ip2_FrictMat_FrictMat_FrictPhys(),
                              
                              #Liquid
                              Ip2_ViscElCapMat_ViscElCapMat_ViscElCapPhys()],
                             
                             
                             [Law2_ScGeom_FrictPhys_CundallStrack(),
                              
                              #liquid
                              Law2_ScGeom_ViscElCapPhys_Basic()]
                             #[Law2_ScGeom_CapillaryPhys_Capillarity],
                             
                             ),
            
            #GlobalStiffnessTimeStepper(),
            NewtonIntegrator(gravity=(0, 0, 0), damping = 0),
            ]


#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
#Cohesion law test
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++


#Init scene -------------------------------------------------------------------

sep_exp = 3 * radius

x_og_pos = 2 * radius
y_exp_pos = 0

fp_id = []
mp_id = []

for cm in range(nbCm) :

    #add particles
    #1
    fp_id.append( O.bodies.append([sphere(center=(0, y_exp_pos, 0),
                                 radius = radius,
                                 material = Models[cm],
                                 fixed = True)])[0] )
    
    
    
    #2
    
    mp_id.append( O.bodies.append([sphere(center=(x_og_pos, y_exp_pos, 0),
                                 radius = radius,
                                 material = Models[cm],
                                 fixed = True)])[0] )
    
    O.bodies[mp_id[cm]].state.vel = [p_speed,0,0]
    
    y_exp_pos += sep_exp
            

#Define Time-step
O.dt = tc_frac * tc


#Energy tracking
O.trackEnergy = True


#Init Var ---------------------------------------------------------------------

ntt_iter = 0
last_time = 0

time = np.full([1],0)

#O.run(1,True)

first = True


#open graphique ctrl
yade.qt.Controller()
qt.View()
    
plot.plots = {'Displacement' : ( Models[0], Models[1], Models[2], Models[3], Models[4], Models[5] )}

# show the plot on the screen, and update while the simulation runs
plot.plot()


###############################################################################
################################### Compute ###################################
###############################################################################

O.run(1000, True)

#save and reload

O.save('./test')

O.load('./test')

print(O.interactions[1,0].geom.penetrationDepth)

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.