If there is anyone who deals with the regrid method after MartinDiehl's (https://github.com/MarDiehl/regrid), Please check my process.

[namhgoo142]

I want to implement heavy deformation simulaton with DAMASK. The regrid method by Martin Diehl is useful. It works well in tow step regriding. However, If I tried to mutiple times, i.e. more than two, the update of orientation in each cell is not working. The third grid's starting orientation is still same as the starting orientation of second step. It should follow the orientation of the last increment of the second step. Here I upload the script file of my process.

###1st process; 1st regriding and two DAMASK_grid runs.

import tempfile
import shutil
import os

import numpy as np
import h5py

import damask

load = 'tensionX'
grid = '20grains16x16x16'
mat = 'material'
grid2 = '20grains16x16x16-2'

cwd = os.getcwd()
#print(wd := tempfile.mkdtemp())
wd = f'{cwd}/wd1'

def new_cells(F_avg,cells):
return (F_avg@cells * np.max(cells/(F_avg@cells))).astype(int)

def regrid_restart(fname_in,fname_out,mapping_flat):
with (h5py.File(fname_in) as f_in, h5py.File(fname_out,'w') as f_out):
f_in.copy('homogenization',f_out)

    f_out.create_group('phase')
    for label in f_in['phase']:
        f_out['phase'].create_group(label)
        F_e0 = np.matmul(f_in['phase'][label]['F'][()],np.linalg.inv(f_in['phase'][label]['F_p'][()]))
        R_e0, V_e0 = damask.mechanics._polar_decomposition(F_e0, ['R','V'])
        f_out['phase'][label].create_dataset('F',data=np.broadcast_to(np.eye(3),(len(mapping_flat),3,3,)))
        f_out['phase'][label].create_dataset('F_e',data=R_e0[mapping_flat])
        f_out['phase'][label].create_dataset('F_p',data=damask.tensor.transpose(R_e0)[mapping_flat])
        f_out['phase'][label].create_dataset('S',data=np.zeros((len(mapping_flat),3,3)))
        for d in f_in['phase'][label]:
            if d in f_out[f'phase/{label}']: continue
            f_out['phase'][label].create_dataset(d,data=f_in['phase'][label][d][()][mapping_flat])

    f_out.create_group('solver')
    for d in ['F','F_lastInc']:
        f_out['solver'].create_dataset(d,data=np.broadcast_to(np.eye(3),np.append(cells_new.prod(),(3,3))))
    for d in ['F_aim', 'F_aim_lastInc']:
        f_out['solver'].create_dataset(d,data=np.eye(3))
    f_out['solver'].create_dataset('F_aimDot',data=np.zeros((3,3)))
    for d in f_in['solver']:
        if d not in f_out['solver']: f_in['solver'].copy(d,f_out['solver'])

'''

docker run --rm --tty --interactive --volume
${PWD}:/wd
--volume ${PWD}/wd1:/output
--env OMP_NUM_THREADS=8
damaskmultiphysics/damask-grid:latest
--load /wd/tensionX.yaml
--geom /wd/20grains16x16x16.vti
--material /wd/material.yaml
-w /output

'''

r = damask.Result(f'{wd}/{grid}{load}{mat}.hdf5')
r.add_IPF_color([0,0,1])
r.export_VTK(target_dir=cwd)
F_avg = np.average(r.view(increments=-1).get('F'),axis=0)

cells_new = new_cells(F_avg,r.cells)
mapping = damask.grid_filters.regrid(r.size,r.view(increments=-1).get('F').reshape(tuple(r.cells)+(3,3)),cells_new)
mapping_flat = mapping.reshape(-1,order='F')

g = damask.GeomGrid.load(f'{grid}.vti')
g.size = [email protected]
g.assemble(mapping).save(f'{wd}/{grid2}.vti')

regrid_restart(f'{wd}/{grid}{load}{mat}restart.hdf5',f'{wd}/{grid2}{load}_{mat}_restart.hdf5',mapping_flat)

r.view(increments=0).export_DADF5(f'{wd}/{grid2}{load}{mat}.hdf5',mapping=mapping)

with h5py.File(f'{wd}/{grid2}{load}{mat}.hdf5','a') as f:
f['geometry'].attrs['size'] = g.size

shutil.copyfile(f'{cwd}/{load}-2.yaml',f'{wd}/{load}.yaml')
shutil.copyfile(f'{wd}/{grid}{load}{mat}.sta',f'{wd}/{grid2}{load}{mat}.sta')

'''
wd1/tensileX.yaml >> f_restart: 300

docker run --rm --tty --interactive --volume
${PWD}:/wd
--volume ${PWD}/wd1:/output
--env OMP_NUM_THREADS=8
damaskmultiphysics/damask-grid:latest
--load /wd/wd1/tensionX.yaml
--geom /wd/wd1/20grains16x16x16-2.vti
--material /wd/material.yaml
-r 140 -w /output

'''

r = damask.Result(f'{wd}/{grid2}{load}{mat}.hdf5').view_less(increments=0)
r.add_IPF_color([0,0,1])
r.export_VTK(target_dir=cwd)

2nd runs

import tempfile
import shutil
import os

import numpy as np
import h5py

import damask

load = 'tensionX'
grid = '20grains16x16x16-2'
mat = 'material'
grid2 = '20grains16x16x16-3'
cwd = os.getcwd()
#print(wd := tempfile.mkdtemp())
wd = f'{cwd}/wd2'

def new_cells(F_avg,cells):
return (F_avg@cells * np.max(cells/(F_avg@cells))).astype(int)

def regrid_restart(fname_in,fname_out,mapping_flat):
with (h5py.File(fname_in) as f_in, h5py.File(fname_out,'w') as f_out):
f_in.copy('homogenization',f_out)

    f_out.create_group('phase')
    for label in f_in['phase']:
        f_out['phase'].create_group(label)
        F_e0 = np.matmul(f_in['phase'][label]['F'][()],np.linalg.inv(f_in['phase'][label]['F_p'][()]))
        R_e0, V_e0 = damask.mechanics._polar_decomposition(F_e0, ['R','V'])
        f_out['phase'][label].create_dataset('F',data=np.broadcast_to(np.eye(3),(len(mapping_flat),3,3,)))
        f_out['phase'][label].create_dataset('F_e',data=R_e0[mapping_flat])
        f_out['phase'][label].create_dataset('F_p',data=damask.tensor.transpose(R_e0)[mapping_flat])
        f_out['phase'][label].create_dataset('S',data=np.zeros((len(mapping_flat),3,3)))
        for d in f_in['phase'][label]:
            if d in f_out[f'phase/{label}']: continue
            f_out['phase'][label].create_dataset(d,data=f_in['phase'][label][d][()][mapping_flat])

    f_out.create_group('solver')
    for d in ['F','F_lastInc']:
        f_out['solver'].create_dataset(d,data=np.broadcast_to(np.eye(3),np.append(cells_new.prod(),(3,3))))
    for d in ['F_aim', 'F_aim_lastInc']:
        f_out['solver'].create_dataset(d,data=np.eye(3))
    f_out['solver'].create_dataset('F_aimDot',data=np.zeros((3,3)))
    for d in f_in['solver']:
        if d not in f_out['solver']: f_in['solver'].copy(d,f_out['solver'])

shutil.copyfile(f'{cwd}/wd1/{grid}{load}{mat}.hdf5',f'{cwd}/wd2/{grid}{load}{mat}.hdf5')
shutil.copyfile(f'{cwd}/wd1/{grid}{load}{mat}restart.hdf5',f'{cwd}/wd2/{grid}{load}{mat}restart.hdf5')
shutil.copyfile(f'{cwd}/wd1/{grid}{load}{mat}.sta',f'{cwd}/wd2/{grid}{load}{mat}.sta')

shutil.copyfile(f'{cwd}/wd1/{grid}.vti',f'{cwd}/wd2/{grid}.vti')

r = damask.Result(f'{cwd}/wd2/{grid}{load}{mat}.hdf5')
#r.add_IPF_color([0,0,1])
#r.export_VTK(target_dir=cwd)
F_avg = np.average(r.view(increments=-1).get('F'),axis=0)

cells_new = new_cells(F_avg,r.cells)
mapping = damask.grid_filters.regrid(r.size,r.view(increments=-1).get('F').reshape(tuple(r.cells)+(3,3)),cells_new)
mapping_flat = mapping.reshape(-1,order='F')

g = damask.GeomGrid.load(f'{wd}/{grid}.vti')
g.size = [email protected]
g.assemble(mapping).save(f'{wd}/{grid2}.vti')

regrid_restart(f'{wd}/{grid}{load}{mat}restart.hdf5',f'{wd}/{grid2}{load}_{mat}_restart.hdf5',mapping_flat)

r.view(increments=0).export_DADF5(f'{wd}/{grid2}{load}{mat}.hdf5',mapping=mapping)

with h5py.File(f'{wd}/{grid2}{load}{mat}.hdf5','a') as f:
f['geometry'].attrs['size'] = g.size

shutil.copyfile(f'{cwd}/{load}-2.yaml',f'{wd}/{load}.yaml')
shutil.copyfile(f'{wd}/{grid}{load}{mat}.sta',f'{wd}/{grid2}{load}{mat}.sta')

'''
N:600
wd2/tensileX.yaml >> f_restart: 600

docker run --rm --tty --interactive --volume
${PWD}:/wd
--volume ${PWD}/wd2:/output
--env OMP_NUM_THREADS=8
damaskmultiphysics/damask-grid:latest
--load /wd/wd2/tensionX.yaml
--geom /wd/wd2/20grains16x16x16-3.vti
--material /wd/material.yaml
-r 340 -w /output

'''

r = damask.Result(f'{wd}/{grid2}{load}{mat}.hdf5').view_less(increments=0)
r.add_IPF_color([0,0,1])
r.export_VTK(target_dir=cwd)