Many errors with PeriodicFlowEngine

Hi, removing particles is not expected with flow engine. I think the simplest solution is to delete the previous FlowEngine after erasing some particles, and replace it with a fresh one.
Probably a good feature would be to let FlowEngine be reset, to avoid such replacement, but the effect is the same.
I hope it helps
Bruno

Thanks Bruno Chareyre, that solved my question.

So how to delete the previous FlowEngine?

Just replace it with a new one in O.engines.
Something like:
O.engines = O.engines[:3]+FlowEngine()+O.engines[4:]

B

Thanks Bruno Chareyre, that solved my question.

Dear Chareyre,
I modified the code in Section B in the way you suggested, but the same error still occurred. Thank you for your attention to this matter, I am honored and grateful for your support and eagerly look forward to your thoughts
This is the revised code in section B：
#B. Activate flow engine
flow.dead=0
flow.defTolerance=-1
flow.meshUpdateInterval=-1
flow.useSolver=0
flow.permeabilityFactor=1
flow.viscosity=10
flow.gradP=Vector3(0,1,0)
flow.duplicateThreshold=0.3
flow.updateTriangulation=True
O.run(2,1)
csdList=flow.getConstrictionsFull()
print(len(O.bodies),len(csdList),'finished')
flow.dead=1
print(O.bodies[10].shape.radius)
for i in range(10,90):
    O.bodies.erase(i)
print(len(O.bodies))
O.engines=O.engines[0:3]+[PeriodicFlowEngine(dead=1,label='flow2')]+O.engines[4:]
O.run(1000,True)
flow2.dead=0
flow2.defTolerance=-1
flow2.meshUpdateInterval=-1
flow2.useSolver=0
flow2.permeabilityFactor=1
flow2.viscosity=10
flow2.gradP=Vector3(0,1,0)
flow2.duplicateThreshold=0.3
flow2.updateTriangulation=True

O.run(2,1)
csdList1=flow2.getConstrictionsFull()
print(len(csdList1))

Do you confirm that the older engine is replaced after this (the indices 3 &4 where just random examples in my post)?
Also, why do you mention deleting particles while your example script does not delete particles?
Do you confirm that no problem occurs when you don't delete particles?
Bruno

Dear Chareyre,
The index of the older engine is 3, so I confirm the older engine is replaced by O.engines=O.engines[0:3]+[PeriodicFlowEngine(dead=1,label='flow2')]+O.engines[4:]
I delete partilces in the code:
for i in range(10,90):
    O.bodies.erase(i)
According to your prompt, I run my code without deleting particles and the following error appears:
CHOLMOD warning: matrix not positive definite. file: ../Supernodal/t_cholmod_super_numeric.c line: 911
something went wrong in Cholesky factorization, use LDLt as fallback this time1
The code is :
##______________ First section, generate sample_________

from __future__ import print_function
from yade import pack, qt, plot
from math import *

nRead=readParamsFromTable(
        ## model parameters
        num_spheres=100,
        targetPorosity= .4,
        confiningPressure=-100000,
        ## material parameters
        compFricDegree=15,#contact friction during the confining phase
        finalFricDegree=30,#contact friction during the deviatoric loading
        young=2e8,
        poisson=.2,
        density=2600,
        alphaKr=7.5,
        alphaKtw=0,
 competaRoll=.22,
        finaletaRoll=.22,
        etaTwist=0,
        normalCohesion=0,
        shearCohesion=0,
        ## fluid parameters
        fluidDensity=1000,
        dynamicViscosity=.001,
        ## control parameters
        damp=0,
        stabilityThreshold=.001,
        ## output specifications
        filename='suffusion',
        unknowOk=True
)

from yade.params.table import *

O.periodic=True
O.cell.hSize=Matrix3(1,0,0, 0,1,0, 0,0,1)
# create materials for spheres
#shear strength is the sum of friction and adhesion, so the momentRotationLaw=True
O.materials.append(CohFrictMat(alphaKr=alphaKr,alphaKtw=alphaKtw,density=density,etaRoll=competaRoll,etaTwist=etaTwist,frictionAngle=radians(compFricDegree),momentRotationLaw=True,normalCohesion=normalCohesion,poisson=poisson,shearCohesion=shearCohesion,young=young,label='spheres'))

# generate particles packing
sp=pack.SpherePack()
sp.makeCloud((0,0,0),(1,1,1),-1,0.3333,num_spheres,False, 0.95,seed=1)
sp.toSimulation(material='spheres')

O.engines=[
        ForceResetter(),
        InsertionSortCollider([Bo1_Sphere_Aabb()]),
        InteractionLoop(
            [Ig2_Sphere_Sphere_ScGeom6D()],
            [Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(label='contact',setCohesionNow=False,setCohesionOnNewContacts=False)],
            [Law2_ScGeom6D_CohFrictPhys_CohesionMoment(useIncrementalForm=True,always_use_moment_law=True)],
 ),
        PeriodicFlowEngine(dead=1,label="flow"),#introduced as a dead engine for the moment, see 2nd section
        GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8),
        PeriTriaxController(label='triax',
            # specify target values and whether they are strains or stresses
            goal=(confiningPressure,confiningPressure,confiningPressure), stressMask=7,
            # type of servo-control, the strain rate isn't determined, it shloud check the unbalanced force
            dynCell=True,maxStrainRate=(10,10,10),
            # wait until the unbalanced force goes below this value
            maxUnbalanced=stabilityThreshold,relStressTol=1e-3,
            doneHook='compactionFinished()'
            ),
        NewtonIntegrator(damping=0)
]

import sys
def compactionFinished():
    # after sample preparation, save the state
    O.save('compactedState'+filename+'.yade.gz')
    print('Compacted state saved', 'porosity', porosity())
    # next time, called python command
    triax.doneHook=''
    O.pause()
O.run()
O.wait()

#B. Activate flow engine
flow.dead=0
flow.defTolerance=-1
flow.meshUpdateInterval=1
flow.useSolver=0
flow.permeabilityFactor=1
flow.viscosity=10
flow.gradP=Vector3(0,1,0)
flow.duplicateThreshold=0.3
flow.updateTriangulation=True
O.run(1000,1)
csdList=flow.getConstrictionsFull()
print(len(O.bodies),len(csdList),'finished')

> I delete partilces in the code:

Indeed, my bad.

> According to your prompt, I run my code without deleting particles and the following error appears

Ok... then that's not about deleting particles.
Since you started from a working script, question then is what you did to make it _not_ work?
Bruno

@Bruno Chareyre
Dear Chareyre,
Thank you for your quick reply! I tried to modify my code, but it failed. I don't know what's wrong with it. I was wondering if you can help me. And I modified the particle size distribution in the FluidCouplingPFV[1] example, then I got the error 'Failed to triangulate body with id=47 Point=-1976.72 1535.93 -26.6935 rad=0.03'. How to eliminate this problem？
Again, I would like to express my most sincere gratitude for the great support!
https://gitlab.com/yade-dev/trunk/blob/master/examples/FluidCouplingPFV/oedometer.py

> I tried to modify my code, but it failed. I don't know what's wrong with it. I was wondering if you can help me.

This really is reminiscent of a sentence in [1] (point 6). :)
So, no, I am not able to help, sorry.
Bruno

[1] https://www.yade-dem.org/wiki/Howtoask