Loading the stored the position of plank, but it wrong

Hello,

please try to provide a MWE [1], where:
- M stands for minimal. Do you have the same problem with one plank? then please provide a script with just one plank..
- W stands for working. I have tried your code with IndentationError: unexpected indent on line 8

> So I storage the position of them, and hope display the rotation process of them, but the reload location is wrong.
> So if reload the location ... the rotation angle does not look right.

store also the rotation (I have not found it in the script)

> the rotation angle does not look right.

as suggested before, use a numerical evidence rather than visual evidence (e.g. print rotation / vertices before loading and after loading)

cheers
Jan

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

To Jan Stránský (honzik) :
Thank you very much. I attach the script as follow. There are two scripts: (1) the simulation to produce the position of 5 windows (2) display the process.

====================================================================================
script (1) Simulation script:

#!/usr/bin/env python
#  encoding: utf-8
from __future__ import print_function
import sys

sys.path.append('~/PycharmProjects/Yade20191229/')
# from yadeImport import *

import math
import numpy as np
import os
from yade.gridpfacet import *
o = Omega()
o.dt = 1e-12
outputDir = 'Output'
checkCounter = 0
numWinds = 5
WinAction = np.zeros((5,1), dtype=int)

## =========== Environment ============##
widthSpace = 6.8 ## width is the distance from the front coal wall to the back
widthHydr = 1.5
lengthSpace = widthHydr * numWinds ## length is the width of 5 windows
highHydr = 3.8
CheckThick = 1
highDummy = 3
colorWind = (0, 1, 1)
angleShield = 50 * 3.1415926 / 180
angleSwingPositive = 15 * 3.1415926 / 180
angleSwingNegtive = 40 * 3.1415926 / 180
lengthShield = 3
lengthTail = 2
windUpperBoundary = 0.9
windLowerBoundary = 0.5019
stateUpperBoundary = highHydr
stateLowerBoundary = windLowerBoundary

positionWind = []
windPosition = np.zeros(5)

shield_y_0 = lengthTail * math.cos(angleShield - angleSwingPositive)
shield_y_1 = shield_y_0 + lengthShield * math.cos(angleShield)
shield_z_0 = highHydr - lengthShield * math.sin(angleShield)
shield_z_1 = highHydr #
wind_y_0 = lengthTail * (math.cos(angleShield - angleSwingPositive) - math.cos(angleShield))
wind_y_1 = shield_y_0
wind_z_0 = highHydr - (lengthShield + lengthTail) * math.sin(angleShield)
wind_z_1 = highHydr - lengthShield * math.sin(angleShield)
HyMat = O.materials.append(FrictMat(frictionAngle=0.1, density=3000, young=2e9))
myGraviaty = -1200
def HydraulicSupport():
    for i in range(0, numWinds):

        temp = [
            Vector3(widthHydr * i, wind_y_0, wind_z_0),
            Vector3(widthHydr * (i + 1), wind_y_0, wind_z_0),
            Vector3(widthHydr * (i + 1), wind_y_1, wind_z_1),
            Vector3(widthHydr * i, wind_y_1, wind_z_1)
        ]
        positionWind.append(temp)

    Wind1 = pack.sweptPolylines2gtsSurface([positionWind[0]], capStart=True, capEnd=True)
    Wind2 = pack.sweptPolylines2gtsSurface([positionWind[1]], capStart=True, capEnd=True)
    Wind3 = pack.sweptPolylines2gtsSurface([positionWind[2]], capStart=True, capEnd=True)
    Wind4 = pack.sweptPolylines2gtsSurface([positionWind[3]], capStart=True, capEnd=True)
    Wind5 = pack.sweptPolylines2gtsSurface([positionWind[4]], capStart=True, capEnd=True)

    IDWind1 = O.bodies.append(pack.gtsSurface2Facets(Wind1, color=colorWind))
    IDWind2 = O.bodies.append(pack.gtsSurface2Facets(Wind2, color=colorWind))
    IDWind3 = O.bodies.append(pack.gtsSurface2Facets(Wind3, color=colorWind))
    IDWind4 = O.bodies.append(pack.gtsSurface2Facets(Wind4, color=colorWind))
    IDWind5 = O.bodies.append(pack.gtsSurface2Facets(Wind5, color=colorWind))

    IDWind = [IDWind1, IDWind2, IDWind3, IDWind4, IDWind5]
    WindList = IDWind1 + IDWind2 + IDWind3 + IDWind4 + IDWind5
    return IDWind, WindList

##---------------------------------------##
def WindowsAction(IDWind):
    global WinAction, windPosition, saveCounter, nEpisode
    RotationW = [RotationW1, RotationW2, RotationW3, RotationW4, RotationW5]

    for nW in range(0, numWinds):
        RotationW[nW].angularVelocity = -0.000001
    SaveProcessLocation(nEpisode, saveCounter)
    saveCounter += 1

def SaveProcessLocation(Episode, iter):
    outputDir = "Output/location/"
    if not os.path.exists(outputDir):
        os.makedirs(outputDir)
    Position = []
    for i in WindList:
        print(o.bodies[i].state.pos)
        temp = [i, o.bodies[i].state.pos]
        Position.append(temp)
    location_name = outputDir + 'wind_' + str(Episode)+ '_' + str(iter)
    np.save(location_name, Position)

def WindowsControl():
    global saveCounter,nEpisode
    nEpisode += 1
    saveCounter = 0

nEpisode = 0
saveCounter = 0
IDWind, WindList = HydraulicSupport()

O.engines = [
    ForceResetter(),
    InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
    InteractionLoop(
        [Ig2_Sphere_Sphere_ScGeom(), Ig2_Facet_Sphere_ScGeom(), Ig2_Wall_Sphere_ScGeom()],
        [Ip2_FrictMat_FrictMat_FrictPhys(),
         Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(xSectionWeibullScaleParameter=0.5,
        xSectionWeibullShapeParameter=0.5,
        weibullCutOffMin=0,
        weibullCutOffMax=10)],
        [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=False, recordMoments=False,label='interactionLaw'),
        Law2_ScGeom_FrictPhys_CundallStrack()]
        ),
    GlobalStiffnessTimeStepper(),

    NewtonIntegrator(gravity=(0, 0, myGraviaty), damping=0.5, label='down'),
    RotationEngine(rotationAxis=(1, 0, 0), rotateAroundZero=True, ids=IDWind[0], zeroPoint=positionWind[0][2],
                   label='RotationW1'),
    RotationEngine(rotationAxis=(1, 0, 0), rotateAroundZero=True, ids=IDWind[1], zeroPoint=positionWind[1][2],
                   label='RotationW2'),
    RotationEngine(rotationAxis=(1, 0, 0), rotateAroundZero=True, ids=IDWind[2], zeroPoint=positionWind[2][2],
                   label='RotationW3'),
    RotationEngine(rotationAxis=(1, 0, 0), rotateAroundZero=True, ids=IDWind[3], zeroPoint=positionWind[3][2],
                   label='RotationW4'),
    RotationEngine(rotationAxis=(1, 0, 0), rotateAroundZero=True, ids=IDWind[4], zeroPoint=positionWind[4][2],
                   label='RotationW5'),

    PyRunner(command="WindowsAction(IDWind)", iterPeriod=200000),
    PyRunner(command="WindowsControl()", iterPeriod=2000000),
]

=====================================================================================
script (2) Display the process

#!/usr/bin/env python
#  encoding: utf-8
from __future__ import print_function
import sys
sys.path.append('~/PycharmProjects/Yade20191229/')
# from yadeImport import *
import math
import numpy as np
import os
from yade.gridpfacet import *
o = Omega()
o.dt = 1e-12
outputDir = 'Output'
checkCounter = 0
numWinds = 5
WinAction = np.zeros((5,1), dtype=int)

## =========== Environment ============##
widthSpace = 6.8 ## width is the distance from the front coal wall to the back
widthHydr = 1.5
lengthSpace = widthHydr * numWinds ## length is the width of 5 windows
highHydr = 3.8
CheckThick = 1
highDummy = 3
colorWind = (0, 1, 1)
angleShield = 50 * 3.1415926 / 180
angleSwingPositive = 15 * 3.1415926 / 180
angleSwingNegtive = 40 * 3.1415926 / 180
lengthShield = 3
lengthTail = 2
windUpperBoundary = 0.9
windLowerBoundary = 0.5019
stateUpperBoundary = highHydr
stateLowerBoundary = windLowerBoundary
positionWind = []
windPosition = np.zeros(5)

shield_y_0 = lengthTail * math.cos(angleShield - angleSwingPositive)
shield_y_1 = shield_y_0 + lengthShield * math.cos(angleShield)
shield_z_0 = highHydr - lengthShield * math.sin(angleShield)
shield_z_1 = highHydr #
wind_y_0 = lengthTail * (math.cos(angleShield - angleSwingPositive) - math.cos(angleShield))
wind_y_1 = shield_y_0
wind_z_0 = highHydr - (lengthShield + lengthTail) * math.sin(angleShield)
wind_z_1 = highHydr - lengthShield * math.sin(angleShield)

HyMat = O.materials.append(FrictMat(frictionAngle=0.1, density=3000, young=2e9))
myGraviaty = -1200
def Relocation():
    global nEpisode, nIter, index
    path = 'Output/location/'
    wind_name = path + 'wind_' + str(nEpisode) + '_' + str(nIter) + '.npy'
    if nIter >= index[nEpisode, 1]:
        nIter = 0
        nEpisode += 1
    else:
        nIter += 1
    if nEpisode >= len(index):
        O.pause()

    locationWind = np.load(wind_name, allow_pickle=True)
    print(locationWind)

    for i in range(0, len(locationWind)):
        n = locationWind[i][0]
        o.bodies[n].state.pos = locationWind[i][1]

def HydraulicSupport():
    for i in range(0, numWinds):
        temp = [
            Vector3(widthHydr * i, wind_y_0, wind_z_0),
            Vector3(widthHydr * (i + 1), wind_y_0, wind_z_0),
            Vector3(widthHydr * (i + 1), wind_y_1, wind_z_1),
            Vector3(widthHydr * i, wind_y_1, wind_z_1)
        ]
        positionWind.append(temp)

    Wind1 = pack.sweptPolylines2gtsSurface([positionWind[0]], capStart=True, capEnd=True)
    Wind2 = pack.sweptPolylines2gtsSurface([positionWind[1]], capStart=True, capEnd=True)
    Wind3 = pack.sweptPolylines2gtsSurface([positionWind[2]], capStart=True, capEnd=True)
    Wind4 = pack.sweptPolylines2gtsSurface([positionWind[3]], capStart=True, capEnd=True)
    Wind5 = pack.sweptPolylines2gtsSurface([positionWind[4]], capStart=True, capEnd=True)

    IDWind1 = O.bodies.append(pack.gtsSurface2Facets(Wind1, color=colorWind))
    IDWind2 = O.bodies.append(pack.gtsSurface2Facets(Wind2, color=colorWind))
    IDWind3 = O.bodies.append(pack.gtsSurface2Facets(Wind3, color=colorWind))
    IDWind4 = O.bodies.append(pack.gtsSurface2Facets(Wind4, color=colorWind))
    IDWind5 = O.bodies.append(pack.gtsSurface2Facets(Wind5, color=colorWind))

    IDWind = [IDWind1, IDWind2, IDWind3, IDWind4, IDWind5]
    WindList = IDWind1 + IDWind2 + IDWind3 + IDWind4 + IDWind5
    return IDWind, WindList

IDWind, WindList = HydraulicSupport()
import os.path
rootdir = "Output/location/"
wind_locations = []
for parent, dirnames, filenames in os.walk(rootdir):
    for filename in filenames:
        if filename[0:4] == 'wind':
            wind_locations.append(filename)
episodeMax = 0

for i in wind_locations:
   n = int(i[5])
   if n >= episodeMax:
       episodeMax = n

index = np.zeros((episodeMax+1, 2))
for i in range(episodeMax+1):
    index[i,0] = i
    for name in wind_locations:
        episode = int(name[5])
        iter = int(name[7])
        if episode == i:
            if (iter) >= index[i, 1]:
                index[i, 1] = iter

nEpisode =0
nIter = 0

O.engines = [
    ForceResetter(),
    InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
    InteractionLoop(
        [Ig2_Sphere_Sphere_ScGeom(), Ig2_Facet_Sphere_ScGeom(), Ig2_Wall_Sphere_ScGeom()],
        [Ip2_FrictMat_FrictMat_FrictPhys(),
         Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(xSectionWeibullScaleParameter=0.5,
        xSectionWeibullShapeParameter=0.5,
        weibullCutOffMin=0,
        weibullCutOffMax=10)],
        [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=False, recordMoments=False,label='interactionLaw'),
        Law2_ScGeom_FrictPhys_CundallStrack()]
        ),
    GlobalStiffnessTimeStepper(),
    NewtonIntegrator(gravity=(0, 0, myGraviaty), damping=0.5, label='down'),
    PyRunner(command="Relocation()", iterPeriod=200000),

]

The position loaded in the second script is not right. That seems the stored location is wrong.
Thank you very much.

Instant replay:
> M stands for minimal. Do you have the same problem with one plank? then please provide a script with just one plank..

please read [1]

A true MWE for this case would be:
- create ONE (or TWO) facet(s)
- move / rotate it a bit, a few iterations is enough
- save it
- load it
does it work ?
- no, then ask your question providing this true MWE
- yes, then add some other minimal part of the code, goto "does it work?"

> The position loaded in the second script is not right

how do you know?
what is "position"?
I have tested your code, and both saved and loaded printed values are the same..

> That seems the stored location is wrong.

seems? either it is or not.. and from the printed values they are correct..

The standard way to view / visualize the simulation is to export data for a dedicated software for visualization, like Paraview.
You can use Yade to visualize the results, of course, but then I don't understand most of the second script, which seems to be "too much simulation" and too similar to the first script to be the visualization script..

cheers
Jan

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

To Jan Stránský:
Thank you very much. I modified the script as follow, there just a plank. In the first It rotates around a axis and the location stored. Then the location load by the second script and display in the 3D show.
The plank consist of two triangles. The 3D show of the first script that the two triangles rotate together, however in the display script, the triangles are separated. As you found the location is the same. So please help me how to modify that.

(1) Rotation script
#!/usr/bin/env python
# encoding: utf-8
from __future__ import print_function
import sys

sys.path.append('~/PycharmProjects/Yade20191229/')
# from yadeImport import *
import numpy as np
import os
from yade.gridpfacet import *
o = Omega()
o.dt = 1e-12
outputDir = 'Output'
positionWind = []

def HydraulicSupport():
    global positionWind
    temp = [Vector3(0, 0, 0),
            Vector3(0, 10, 0),
            Vector3(10, 10, 0),
            Vector3(10, 0, 0)]
    positionWind.append(temp)
    Wind1 = pack.sweptPolylines2gtsSurface([positionWind[0]], capStart=True, capEnd=True)
    IDWind = O.bodies.append(pack.gtsSurface2Facets(Wind1))
    return IDWind

def Ground():
    O.bodies.append(utils.wall(position=-10, sense=0, axis=2, material=-1))

##---------------------------------------##
def WindowsAction(IDWind):
    global WinAction, windPosition, saveCounter, nEpisode,RotationW1
    RotationW1.angularVelocity = -0.000001
    SaveProcessLocation(nEpisode, saveCounter)
    saveCounter += 1

def SaveProcessLocation(Episode, iter):
    outputDir = "Output/location/"
    if not os.path.exists(outputDir):
        os.makedirs(outputDir)
    Position = []

    for i in IDWind:
        print(o.bodies[i].state.pos)
        temp = [i, o.bodies[i].state.pos]
        Position.append(temp)
    location_name = outputDir + 'wind_' + str(Episode)+ '_' + str(iter)
    np.save(location_name, Position)

def WindowsControl():
    global saveCounter, nEpisode
    nEpisode += 1
    saveCounter = 0

nEpisode = 0
saveCounter = 0
IDWind = HydraulicSupport()
Ground()

O.engines = [
    ForceResetter(),
    InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
    InteractionLoop(
        [Ig2_Sphere_Sphere_ScGeom(), Ig2_Facet_Sphere_ScGeom(), Ig2_Wall_Sphere_ScGeom()],
        [Ip2_FrictMat_FrictMat_FrictPhys(),
         Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(xSectionWeibullScaleParameter=0.5,
        xSectionWeibullShapeParameter=0.5,
        weibullCutOffMin=0,
        weibullCutOffMax=10)],
        [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=False, recordMoments=False,label='interactionLaw'),
        Law2_ScGeom_FrictPhys_CundallStrack()]
        ),
    GlobalStiffnessTimeStepper(),
    NewtonIntegrator(gravity=(0, 0, 9.8), damping=0.5, label='down'),
    RotationEngine(rotationAxis=(1, 0, 0), rotateAroundZero=True, ids=IDWind, zeroPoint=positionWind[0][0],
                   label='RotationW1'),
    PyRunner(command="WindowsAction(IDWind)", iterPeriod=200000),
    PyRunner(command="WindowsControl()", iterPeriod=2000000),
]
# o.run(1200000000)

================================================================
(2) Display script
from __future__ import print_function
import sys
sys.path.append('~/PycharmProjects/Yade20191229/')
# from yadeImport import *
import numpy as np
import os
from yade.gridpfacet import *
o = Omega()
o.dt = 1e-12

positionWind=[]

def HydraulicSupport():
    global positionWind
    temp = [Vector3(0, 0, 0),
            Vector3(0, 10, 0),
            Vector3(10,10, 0),
            Vector3(10,0, 0)]
    positionWind.append(temp)
    Wind1 = pack.sweptPolylines2gtsSurface([positionWind[0]], capStart=True, capEnd=True)
    IDWind = O.bodies.append(pack.gtsSurface2Facets(Wind1))
    return IDWind

def Ground():
    O.bodies.append(utils.wall(position=-10, sense=0, axis=2, material=-1))

def Relocation():
    global nEpisode, nIter, index
    path = 'Output/location/'
    wind_name = path + 'wind_' + str(nEpisode) + '_' + str(nIter) + '.npy'

    if nIter >= index[nEpisode, 1]:
        nIter = 0
        nEpisode += 1
    else:
        nIter += 1
    if nEpisode >= len(index):
        O.pause()

    locationWind = np.load(wind_name, allow_pickle=True)
    print(nEpisode,nIter)
    print(locationWind)
    print('--------------')
    for i in range(0, len(locationWind)):
        n = locationWind[i][0]
        o.bodies[n].state.pos = locationWind[i][1]

import os.path
rootdir = "Output/location/"
wind_locations = []

for parent, dirnames, filenames in os.walk(rootdir):
    for filename in filenames:
        if filename[0:4] == 'wind':
            wind_locations.append(filename)

episodeMax = 0
endNumber = 100000000
for i in wind_locations:
   n = int(i[5])
   if n >= episodeMax:
       episodeMax = n

index = np.zeros((episodeMax+1, 2))
for i in range(episodeMax+1):
    index[i,0] = i
    for name in wind_locations:
        episode = int(name[5])
        iter = int(name[7])
        if episode == i:
            if (iter) >= index[i, 1]:
                index[i, 1] = iter

IDWind = HydraulicSupport()
Ground()
nEpisode =0
nIter = 0

O.engines = [
    ForceResetter(),
    InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
    InteractionLoop(
        [Ig2_Sphere_Sphere_ScGeom(), Ig2_Facet_Sphere_ScGeom(), Ig2_Wall_Sphere_ScGeom()],
        [Ip2_FrictMat_FrictMat_FrictPhys(),
         Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(xSectionWeibullScaleParameter=0.5,
        xSectionWeibullShapeParameter=0.5,
        weibullCutOffMin=0,
        weibullCutOffMax=10)],
        [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=False, recordMoments=False,label='interactionLaw'),
        Law2_ScGeom_FrictPhys_CundallStrack()]
        ),
    GlobalStiffnessTimeStepper(),

    NewtonIntegrator(gravity=(0, 0, 9.8), damping=0.5, label='down'),
    PyRunner(command="Relocation()", iterPeriod=100000),
]
# o.run(1200000000, True)

Thank you very much

Sorry, but this is still too far from MWE. If you want help from us, you have to do some work on your side first.

If your problem is saving/loading a facet, then the MWE should be about saving and loading one facet:

### 1st script
f = facet(((0,0,0),(1,0,0),(0,1,0)))
O.bodies.append(f)
f.state.vel = (.1,0,0)
f.state.angVel = (0,0,.1)
O.run(10,True)
numpy.save("f",[[f.id,f.state.pos]])
###

### 2nd script
f = facet(((0,0,0),(1,0,0),(0,1,0)))
O.bodies.append(f)
[[i,pos]] = numpy.load("f.npy",allow_pickle=True)
f.state.pos = pos
###

is this the same problem you are facing?

cheers
Jan

To Jan Stránský :
Thank you very much. I am sorry that the script can not meet MWE. My script is similar with yours. The facet rotate and restore the location , then display the rotation process by reload the location. How ever my facet is a quadrangle, and consist with two triangles.
      I used the " pack.sweptPolylines2gtsSurface([positionWind[0]], capStart=True, capEnd=True)" to create the plank, but I do not know why the quadrangle consists with two triangles. That is source of my problem. In the rotation process, the quadrangle moving with the two triangles sticks together in the 3D show. However, If I read the location from the storage, in the 3D show, the two triangles separated.

    I desire your help. so I modified the script as follow, I hope it can meet the MWE

(1) Rotation process:

positionWind = []
temp = [Vector3(0, 0, 0),
        Vector3(0, 10, 0),
        Vector3(10, 10, 0),
        Vector3(10, 0, 0)]
positionWind.append(temp)
Wind1 = pack.sweptPolylines2gtsSurface([positionWind[0]], capStart=True, capEnd=True)
IDWind = O.bodies.append(pack.gtsSurface2Facets(Wind1))
saveCounter = 0
##---------------------------------------##
def WindowsAction_SaveLoaction(IDWind):
    global WinAction, windPosition, saveCounter, RotationW1
    RotationW1.angularVelocity = -0.000001
    outputDir = "Output/location/"
    if not os.path.exists(outputDir):
        os.makedirs(outputDir)
    Position = []
    for i in IDWind:
        temp = [i, o.bodies[i].state.pos]
        Position.append(temp)
    location_name = outputDir + 'wind' + '_' + str(saveCounter)
    np.save(location_name, Position)
    saveCounter += 1

O.engines = [
    ForceResetter(),
    GlobalStiffnessTimeStepper(),
    NewtonIntegrator(gravity=(0, 0, 9.8), damping=0.5, label='down'),
    RotationEngine(rotationAxis=(1, 0, 0), rotateAroundZero=True, ids=IDWind, zeroPoint=positionWind[0][0],
                   label='RotationW1'),
    PyRunner(command="WindowsAction_SaveLoaction(IDWind)", iterPeriod=200000),
]

==========================================
(2) Read the location

positionWind = []
temp = [Vector3(0, 0, 0),
        Vector3(0, 10, 0),
        Vector3(10, 10, 0),
        Vector3(10, 0, 0)]
positionWind.append(temp)
Wind1 = pack.sweptPolylines2gtsSurface([positionWind[0]], capStart=True, capEnd=True)
IDWind = O.bodies.append(pack.gtsSurface2Facets(Wind1))

saveCounter = 0
import os.path
rootdir = "Output/location/"
wind_locations = []
for parent, dirnames, filenames in os.walk(rootdir):
    for filename in filenames:
        if filename[0:4] == 'wind':
            wind_locations.append(filename)
fileNum = len(filenames)
saveCounter=0

def Display():
    global fileNum,saveCounter
    path = 'Output/location/'
    wind_name = path + 'wind_' + str(saveCounter) + '.npy'
    locationWind = np.load(wind_name, allow_pickle=True)
    print(locationWind)
    print('--------------')
    for i in range(0, len(locationWind)):
        n = locationWind[i][0]
        o.bodies[n].state.pos = locationWind[i][1]
    saveCounter +=1

O.engines = [
    ForceResetter(),
    GlobalStiffnessTimeStepper(),
    NewtonIntegrator(gravity=(0, 0, 9.8), damping=0.5, label='down'),
    PyRunner(command="Display()", iterPeriod=200000),
]

Thank you very much

> my facet is a quadrangle, and consist with two triangles.
> I used the " pack.sweptPolylines2gtsSurface([positionWind[0]], capStart=True, capEnd=True)" to create the plank, but I do not know why the quadrangle consists with two triangles. That is source of my problem.

sweptPolylines2gtsSurface - gts means GNU **Triangulated** Surface [1], therefore the result are triangles
Yade facet = triangle. Yade currently does not support more general planar shapes.
You can use box [2] with zero or "small" thick if your aim is a rectangle.

> My script is similar with yours.

Ok, now we are targeting the source of the problem :-)
Yade Body has a State, which has position (pos), but also an orientation (ori). Depending on the shape, to be able to restore the position AND orientation, you should save state.pos, but also state.ori.

> I modified the script as follow, I hope it can meet the MWE

getting closer, thanks :-) but I will update my MWE:
### script 1
f = facet(((0,0,0),(1,0,0),(0,1,0)))
O.bodies.append(f)
f.state.vel = (.1,0,0)
f.state.angVel = (0,0,.1)
O.run(10,True)
data = []
for b in O.bodies:
   axis,angle = b.state.ori.toAxisAngle()
   data.append([b.id,b.state.pos,axis,angle])
numpy.save("f",data)
### script 2
f = facet(((0,0,0),(1,0,0),(0,1,0)))
O.bodies.append(f)
print(numpy.load("f.npy",allow_pickle=True))
[[i,pos,axis,angle]] = numpy.load("f.npy",allow_pickle=True)
f.state.pos = pos
f.state.ori = Quaternion(axis,angle)
###

a few notes:
- as pointed above, consider an external program for visualization
- instead of pos or ori, you can save vertices of facets, which is somehow more general than just pos and ori and depending on the geometry to be hard-coded (similarly to other shapes e.g. you switch to box)

cheers
Jan

[1] http://gts.sourceforge.net/
[2] https://yade-dem.org/doc/yade.utils.html#yade.utils.box