Difference between revisions of "Contrib:KeesWouters/Homard/lshape/loop"
Keeswouters (Talk | contribs) m (→'''The geometry and the mesh of the construction''') |
Keeswouters (Talk | contribs) |
||
Line 27: | Line 27: | ||
* define number of refinement iterations (whether or not fixed) | * define number of refinement iterations (whether or not fixed) | ||
− | * | + | * for counter in [0,1,2,3]: |
* ...define model | * ...define model | ||
* ...apply material to model | * ...apply material to model | ||
Line 41: | Line 41: | ||
The refinement command itself is as follows: | The refinement command itself is as follows: | ||
− | + | ...... | |
− | + | counterp1 = counter+1; | |
− | + | meshcp1=CO('mesh_%d' %counterp1) | |
− | + | MACR_ADAP_MAIL(MAILLAGE_N=mesh[counter], | |
− | + | MAILLAGE_NP1=meshcp1, | |
− | + | ADAPTATION='RAFFINEMENT', | |
− | + | RESULTAT_N=result[counter], | |
− | + | INDICATEUR='ERRE_ELEM_SIGM', | |
− | + | NOM_CMP_INDICA='ERREST', | |
− | + | CRIT_RAFF_PE=0.20, | |
− | + | QUALITE='OUI', | |
− | + | CONNEXITE='OUI', | |
+ | TAILLE='OUI',); | ||
− | The command takes the current mesh '' | + | mesh[counter+1] = meshcp1; |
+ | |||
+ | |||
+ | The command takes the current mesh ''mesh[counter]'' as input '''MAILLAGE_N=curmesh''' and the output mesh is '''MAILLAGE_NP1=meshcp1'''. Again, the meshcp1 is not defined (in the Python sense) so we need to use the CO('..') operator on forehand: | ||
+ | counterp1 = counter+1; | ||
+ | meshcp1=CO('mesh_%d' %counterp1) | ||
+ | |||
+ | Now suppose counter=1, counterp1=2 then meshcp1 then evaluates to '''CO('mesh_2') = mesh_2'''. This corresponds to mesh[2]. | ||
+ | |||
+ | Further keywords are the ADAPTATION='RAFFINEMENT' together with CRIT_RAFF_PE=0.10. In this case the 0.10 or 10 % of the elements with the largest error estimate (based on stress error indicator INDICATEUR='ERRE_ELEM_SIGM', NOM_CMP_INDICA='ERREST') are refined. The other elements remain untouched. | ||
+ | At the end we need to assign the next element of array mesh to the newly defined meshcp1: mesh[counter+1] = meshcp1. | ||
==='''The detailed command of the refinement'''=== | ==='''The detailed command of the refinement'''=== |
Revision as of 00:12, 23 February 2010
Contents
Mesh refinement using Python loop facility
[under construction - just started 2010-02-22 ... ]
This contribution is based on the mesh refinement described here. Two mesh refinements based on three calculations are given there. Since the three calculations are basically the same, it is convenient to use the Python looping facility. So now we will perform four calculations and three mesh refinements. The iteration is controlled by the a Python for loop.
The geometry
The geometry and the mesh of the construction
The mesh refinement is based on Homard. As usual the geometry is quite simple and is defined in the link given above.
Again, the mesh consists of coarse, tetrehedral elements to start with.
The load of the construction
The L shape is fixed in all directions at the bottom area Afix. On area at opposite end of L is denoted Apres. For simplicity, on this area a displacement of 0.1 mm is given in the y direction.
The general idea of the refinement
The general idea of the Homard refinement is to start with a rather course mesh and update the mesh according to a criterium. In this case we use the stress error estimate ERRE_ELEM_SIGM and its component ERREST. For more details on this error estimation see [....]. So to start with, a standard load case will be calculated and an error estimation on all the elements will be carried out. Again depending on a criterium the mesh will be adapted. In this example a certain percentage of the elements with the highest error estimates will be refined. A new calculation with the adapted mesh can take place. For this calculation all the models, loads and parameters need to be rebuild. So it is most suitable to define geometrical entities (GROUP_MA) in the mesh in stead of nodes.
To summarise:
- initial some parameters
- read initial mesh
- define material properties
- define number of refinement iterations (whether or not fixed)
- for counter in [0,1,2,3]:
- ...define model
- ...apply material to model
- ...apply load to model
- ...perform calculation
- ...determine element and node parameters
- ...(write output data)
- ...refinement of mesh MACR_ADAP_MAIL
The error estimation is determined in the CALC_ELEM command, by the option 'ERRE_ELEM_SIGM':
result...=CALC_ELEM(....,OPTION=(....,'ERRE_ELEM_SIGM',),);
The refinement command itself is as follows:
...... counterp1 = counter+1; meshcp1=CO('mesh_%d' %counterp1) MACR_ADAP_MAIL(MAILLAGE_N=mesh[counter], MAILLAGE_NP1=meshcp1, ADAPTATION='RAFFINEMENT', RESULTAT_N=result[counter], INDICATEUR='ERRE_ELEM_SIGM', NOM_CMP_INDICA='ERREST', CRIT_RAFF_PE=0.20, QUALITE='OUI', CONNEXITE='OUI', TAILLE='OUI',);
mesh[counter+1] = meshcp1;
The command takes the current mesh mesh[counter] as input MAILLAGE_N=curmesh and the output mesh is MAILLAGE_NP1=meshcp1. Again, the meshcp1 is not defined (in the Python sense) so we need to use the CO('..') operator on forehand:
counterp1 = counter+1; meshcp1=CO('mesh_%d' %counterp1)
Now suppose counter=1, counterp1=2 then meshcp1 then evaluates to CO('mesh_2') = mesh_2. This corresponds to mesh[2].
Further keywords are the ADAPTATION='RAFFINEMENT' together with CRIT_RAFF_PE=0.10. In this case the 0.10 or 10 % of the elements with the largest error estimate (based on stress error indicator INDICATEUR='ERRE_ELEM_SIGM', NOM_CMP_INDICA='ERREST') are refined. The other elements remain untouched. At the end we need to assign the next element of array mesh to the newly defined meshcp1: mesh[counter+1] = meshcp1.
The detailed command of the refinement
In this case we perform two refinements of the mesh, ie. three standard calculations are performed. Each refinement of the mesh is based on a parameter of the previous calculation.
The standard calculation
This set of commands is a standard calculation. The only command that should be added for a complete calculation is FIN():
DEBUT();
#Read MED mesh file mesh0=LIRE_MAILLAGE(UNITE=20,FORMAT='MED',NOM_MED='mesh0',INFO_MED=2,INFO=1,);
#Assign to the Model model0=AFFE_MODELE(MAILLAGE=mesh0,AFFE=_F(TOUT='OUI', PHENOMENE='MECANIQUE',MODELISATION='3D',),);
#Define material properties steel=DEFI_MATERIAU(ELAS=_F(E=2.1e5,NU=.28,),);
#Assign material to mesh Amat0=AFFE_MATERIAU(MAILLAGE=mesh0,AFFE=_F(TOUT='OUI',MATER=steel,),);
# define boundary conditions and loads Load0=AFFE_CHAR_MECA(MODELE=model0, DDL_IMPO=(_F(GROUP_MA='Afix',DX=0.0,DY=0.0,DZ=0.0,), _F(GROUP_MA='Apres',DX=0.0,DY=0.1,DZ=0.0,),),);
result0=MECA_STATIQUE(MODELE=model0,CHAM_MATER=Amat0, EXCIT=_F(CHARGE=Load0,),),);
result0=CALC_ELEM(reuse =result0,MODELE=model0,RESULTAT=result0, TOUT='OUI',TYPE_OPTION='TOUTES', OPTION=('EQUI_ELNO_SIGM','SIEF_ELNO_ELGA','ERRE_ELEM_SIGM',),);
result0=CALC_NO(reuse =result0,RESULTAT=result0,OPTION=('EQUI_NOEU_SIGM',),);
IMPR_RESU(FORMAT='MED',UNITE=80,RESU=_F(MAILLAGE=mesh0,RESULTAT=result0, NOM_CHAM=('EQUI_ELNO_SIGM','EQUI_NOEU_SIGM','ERRE_ELEM_SIGM','DEPL',),),);
The first refinement - on mesh0
The mesh0 is refined. The refined mesh is available under the name mesh1. Now exactly the same, standard calculation is being performed. The models, loads and material properties first need to be based on mesh1. So also the assingments for this new mesh need to be carried out. The results are now written to unit 81 in stead of unit 80.
MACR_ADAP_MAIL(MAILLAGE_N=mesh0, MAILLAGE_NP1=CO('mesh1'), ADAPTATION='RAFFINEMENT', RESULTAT_N=result0, INDICATEUR='ERRE_ELEM_SIGM', NOM_CMP_INDICA='ERREST', CRIT_RAFF_PE=0.30, QUALITE='OUI', CONNEXITE='OUI', TAILLE='OUI',);
model1=AFFE_MODELE(MAILLAGE=mesh1,AFFE=_F(TOUT='OUI', PHENOMENE='MECANIQUE',MODELISATION='3D',),);
Load1=AFFE_CHAR_MECA(MODELE=model1, DDL_IMPO=(_F(GROUP_MA='Afix',DX=0.0,DY=0.0,DZ=0.0,), _F(GROUP_MA='Apres',DX=0.0,DY=0.1,DZ=0.0,),),);
#Assign material to mesh Amat1=AFFE_MATERIAU(MAILLAGE=mesh1,AFFE=_F(TOUT='OUI',MATER=steel,),);
result1=MECA_STATIQUE(MODELE=model1,CHAM_MATER=Amat1, EXCIT=(_F(CHARGE=Load1,),),);
result1=CALC_ELEM(reuse =result1,MODELE=model1,RESULTAT=result1, TOUT='OUI',TYPE_OPTION='TOUTES', OPTION=('EQUI_ELNO_SIGM','SIEF_ELNO_ELGA','ERRE_ELEM_SIGM',),);
result1=CALC_NO(reuse =result1,RESULTAT=result1,OPTION=('EQUI_NOEU_SIGM',),);
IMPR_RESU(FORMAT='MED',UNITE=81, RESU=_F(MAILLAGE=mesh1,RESULTAT=result1, NOM_CHAM=('EQUI_ELNO_SIGM','SIGM_NOEU_DEPL','DEPL','SIEF_ELNO_ELGA','ERRE_ELEM_SIGM',),),);
The second refinement - on mesh1
Same procedure as before but now carried out on mesh1 and the output is the refined mesh2. Of course, the models, loads and material properties need to be based on mesh2 now. The results are now written to unit 82.
MACR_ADAP_MAIL(MAILLAGE_N=mesh1, MAILLAGE_NP1=CO('mesh2'), ADAPTATION='RAFFINEMENT', RESULTAT_N=result1, INDICATEUR='ERRE_ELEM_SIGM', NOM_CMP_INDICA='ERREST', CRIT_RAFF_PE=0.30, QUALITE='OUI', CONNEXITE='OUI', TAILLE='OUI',);
model2=AFFE_MODELE(MAILLAGE=mesh2,AFFE=_F(TOUT='OUI', PHENOMENE='MECANIQUE',MODELISATION='3D',),); Load2=AFFE_CHAR_MECA(MODELE=model2, DDL_IMPO=(_F(GROUP_MA='Afix',DX=0.0,DY=0.0,DZ=0.0,), _F(GROUP_MA='Apres',DX=0.0,DY=0.1,DZ=0.0,),),);
#Assign material to mesh Amat2=AFFE_MATERIAU(MAILLAGE=mesh2,AFFE=_F(TOUT='OUI',MATER=steel,),);
result2=MECA_STATIQUE(MODELE=model2,CHAM_MATER=Amat2, EXCIT=(_F(CHARGE=Load2,),),);
result2=CALC_ELEM(reuse =result2,MODELE=model2,RESULTAT=result2, TOUT='OUI',TYPE_OPTION='TOUTES', OPTION=('EQUI_ELNO_SIGM','SIEF_ELNO_ELGA','ERRE_ELEM_SIGM',),);
result2=CALC_NO(reuse =result2,RESULTAT=result2,OPTION=('EQUI_NOEU_SIGM',),);
IMPR_RESU(FORMAT='MED',UNITE=82, RESU=_F(MAILLAGE=mesh2,RESULTAT=result2, NOM_CHAM=('EQUI_ELNO_SIGM','SIGM_NOEU_DEPL','DEPL','SIEF_ELNO_ELGA','ERRE_ELEM_SIGM',),),);
FIN(FORMAT_HDF='OUI',);
The results of the refinement
The numbers of nodes and quadratic elements (tetrahedrons) of the meshes:
30% 70% 30%-->new factor nodes elements changed unchanged elements change 4224 2064 620 1444 16346 9636 2890 6746 8192 13 (=8192 / 620) 65845 42537 12760 29777 35791 12 (=35791/2890)
So, from a single tetrahedron that is refined, 12 to 13 tetrahedrons are generated for this mesh.
In the graph below, the original mesh (top left) and the two successive refinements based on the discussed criteria are depicted (bottom right last iteration).
The picture below shows the stress error estimation 'ERREST'. The first refinement is based on this field: mesh0 --> mesh1.
The error estimation field 'ERROR-ELEM_SIGMA' containes 10 components:
ERREST, NUEST, SIGCAL, TERMRE,TERMR2,TERMRENO,TERMN2,TERMSA,TERMS2, TAILLE.
Any of these components can be used in the component indicator, eg. NOM_CMP_INDICA='ERREST'
The ASTK input files
The picture below shows the input and output files defined in ASTK.
Download files
This zip contains the following files
- python geometry and mesh files (load in Salome by File --> Load script (cntrl T)) and right select refresh (F5) in the object browser). Export the med file in the mesh module under mesh0.med for further processing by Code-Aster, controlled by ASTK.
- command files for Code-Aster
- astk file for file control by ASTK
- export file, generated by ASTK.
The result files resi0.med, resi1.med and resi2.med, can be viewed in the post processor module of Salome by File --> Import --> resix.med or cntrl I in the object browser.
Links
More information about the mesh refinement by Homard can be find here:
pdf file: http://www.code-aster.org/V2/doc/default/man_u/u7/u7.03.01.pdf
Homard website: http://www.code-aster.org/outils/homard
That's it for now.
Kees