search

thelfer / MFrontGallery

This project shows how to create a compilation project for MFront properties, behaviours and models
https://thelfer.github.io/MFrontGallery/web/index.html
20 stars 13 forks source link

Wood material repo #2

Open lore-rip opened 2 years ago

lore-rip commented 2 years ago

I would like to ask to include behaviors for wood in the repository. Here attached the first one that is an orthotropic Maxwell viscoelastic law, described in the attached paper. Thanks Lorenzo ortho_Maxwell.zip 3D_visco-elastic_modeling_Froidevaux_al.pdf

thelfer commented 2 years ago

Hi Lorenzo, Thanks for this contribution. The model is well written. It would be nice if we could have unit non-regression tests associated with it. Could you provide some ? Regards, Thomas

lore-rip commented 2 years ago

Hi Thomas, here a regression test featuring a creep test for the three anatomical directions of wood, Longitudinal, Radial, Tangential. An unitary load (1MPa) is applied in each direction. The material law coefficient passed through code_aster are: MFRONT=_F(LISTE_COEF=(10285, 369, 635, 0.448, 0.292, 0.35, 536, 571, 25, 329, 220, 384, 52, 37, 22, 66, 75, 97, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 66, 57, 57, 40, 36, 35, 3, 2, 2, 100, 1000, 10000))) They are representing Poplar wood.

The numerical solution is made by code_aster 14.6 with the command SIMU_POINT_MAT

Thank you

Lorenzo radial_creep.txt long_creep.txt tang_creep.txt

thelfer commented 2 years ago

Hi Lorenzo, Here is a first comparison between MTest and code_aster. mtest_aster. I performed the tests in small strain. Thomas

thelfer commented 2 years ago

Hi Lorenzo,

Created this branch to work on this issue: https://github.com/thelfer/MFrontGallery/tree/th/work_issue_2

I added a first implementation in small strain called OrthotropicGeneralizedMaxwell in the generic-behaviours/viscoelasticity directory as well as an associated test.

I then created a second implementation called PoplarOrthotropicGeneralizedMaxwell_2021 dedicated to poplar in the materials/Wood/behaviours directory. All the material properties have been turned into parameters.

The common part of both implementations is in a file called OrthotropicGeneralizedMaxwell-core.mfront in generic-behaviours/viscoelasticity

There are two things to discuss now:

lore-rip commented 2 years ago

Hi Thomas, Thank you very much. Following the pipe convention the coefficient for wood would be MFRONT=_F(LISTE_COEF=(635, 10285, 369, 0.029, 0.42, 0.165, 786, 838, 114, 66, 75, 97, 329, 220, 384, 37, 22, 66, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 66, 57, 57, 40, 36, 35, 3, 2, 2, 100, 1000, 10000))) The values are literature based so poplar. 2021 is perfect

lore-rip commented 2 years ago

And Here the corrected unit test with SigmaXX=-1Mpa, in radial direction

I would ask you to correct the date in the material behavior, I left a wrong one in the text Many thanks

radial_creep.txt

thelfer commented 2 years ago

Thank you very much. Following the pipe convention the coefficient for wood would be MFRONT=_F(LISTE_COEF=(635, 10285, 369, 0.029, 0.42, 0.165, 786, 838, 114...

Why aren't the shear modulus just a permutation of the previous values ?

lore-rip commented 2 years ago

Because I made an error in the first list. I apologize.

thelfer commented 2 years ago

Hi Lorenzo, Unfortunately, I can't reproduce your reference results with this new parameters. Would you try to run the attached test ? Thomas

$ mfront --obuild --interface=aster OrthotropicGeneralizedMaxwell.mfront 
$ mtest OrthotropicGeneralizedMaxwellRadialTest.mtest

RP.zip

thelfer commented 2 years ago

Hi @lore-rip, My bad, I kept the old reference files. Here is the same test. Results are ok if I switch EYY and EZZ compared to SIMU_POINT_MAT. Is it normal ?

RP.zip

lore-rip commented 2 years ago

Hi Thomas, yes it is normal and correct. Thank you

thelfer commented 2 years ago

Hi @lore-rip Perfect. Would you try to adapt the two other tests ? Thanks, Thomas

lore-rip commented 2 years ago

Hi Thomas, here attached the two tests thank you tang_creep.txt long_creep.txt

thelfer commented 2 years ago

Hi, Both tests are now in the th/work_issue_2 branch. We shall now work on the documentation of the behaviour. To me, it barely boils down reproducing paragraph 2.2 of Froidevaux' paper. I'll make you a first proposal. Best, Thomas

lore-rip commented 2 years ago

Hi, Here attached the .mfront file made for the cylindrical elasticity. The file works here stand alone, where the angle phi is a material property, as well as the moisture content.

Three regression tests are included with material rotation of 0°,45° and 90° and unitary moisture content variation The tests as usual are done with code_aster SIMU_POINT_MAT The three input are MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 114., 0.00064, 0.00013, 0.00115, 1., 0))) MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 114., 0.00064, 0.00013, 0.00115, 1., pi/4))) MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 114., 0.00064, 0.00013, 0.00115, 1., pi/2))

All the best

Lorenzo

elas_cyl.zip )

thelfer commented 2 years ago

@lore-rip. Just one remark, why did you declare the moisture as a material property and not an external state variable.

thelfer commented 2 years ago

@lore-rip Sorry, but I can't reproduce the results of your test. Please find attached a test for phi=0 which uses your implementation. Would you have a look at it ?

elas_cyl.zip

LoreRiparbelli commented 2 years ago

Just one remark, why did you declare the moisture as a material property and not an external state variable.

Hi Yes it is the correct way indeed but it had problems with SIMU_POINT_MAT in code aster Both moisture and phi are external state variable generated a priori of the mechanical analysis

lore-rip commented 2 years ago

@lore-rip Sorry, but I can't reproduce the results of your test. Please find attached a test for phi=0 which uses your implementation. Would you have a look at it ?

elas_cyl.zip

Hi, I noticed that some values from mtest were not correct, I changed my input to match your mtest. Could you check again? cyl_test1.resu.zip

thelfer commented 2 years ago

@lore-rip It works now. What did you change ?

LoreRiparbelli commented 2 years ago

@thelfer I changed the imposed sigma value, the last shear modulus and added the 10 steps. I think that the first one is the only affecting the result. Please read the private mail I sent you. Many thanks

thelfer commented 2 years ago

@lore-rip Indeed, 10 steps shall not change anything.

For the imposed sigma value, I just read it in the first results that you sent me...

For the shear modulus, I used the values that you gave me:

MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 114., 0.00064, 0.00013, 0.00115, 1., 0))) MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 114., 0.00064, 0.00013, 0.00115, 1., pi/4))) MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 114., 0.00064, 0.00013, 0.00115, 1., pi/2))

I did not (yet) receive your private mail, or are you taking about a previous one ?

LoreRiparbelli commented 2 years ago

For the shear modulus, I used the values that you gave me:

@thelfer as far as I understand you used a Shear Modulus13 of 113 instead of 114

I sent the mail again, please confirm to me that if it is not arriving, as it is quite important

thelfer commented 2 years ago

You're right. I got your mail.

LoreRiparbelli commented 2 years ago

Thank you very much!! :)

thelfer commented 2 years ago

@lore-rip Would you have a look at this directory: https://github.com/thelfer/MFrontGallery/tree/master/generic-behaviours/elasticity. Could you send me the results of the tests for angles pi/4 and pi/2 ?

lore-rip commented 2 years ago

@thelfer running your updated script I obtain this error in the command of code_aster (15 and 14) that build the library:

CREA_LIB_MFRONT(DEBUG='NON', UNITE_LIBRAIRIE=7, UNITE_MFRONT=5)

Treating target : all In file included from asterWoodOrthotropicElasticity_2022.cxx:11:0: fort.5:43:1: error: ‘tmatrix’ does not name a type compilation terminated due to -Wfatal-errors. Makefile.mfront:37: recipe for target 'asterWoodOrthotropicElasticity_2022.o' failed make: *** [asterWoodOrthotropicElasticity_2022.o] Error 1 callMake: can't build target 'all' libraries building went wrong

thelfer commented 2 years ago

@lore-rip I updated to make it work with earlier versions of TFEL.

lore-rip commented 2 years ago

@thelfer it seems that a problem with the variable dw is arising

CREA_LIB_MFRONT(DEBUG='NON', UNITE_LIBRAIRIE=7, UNITE_MFRONT=5)

Treating target : all In file included from asterWoodOrthotropicElasticity_2022.cxx:11:0: fort.5: In member function ‘tfel::material::WoodOrthotropicElasticity_2022<(tfel::material::ModellingHypothesis::Hypothesis)6u, Type, false>::IntegrationResult tfel::material::WoodOrthotropicElasticity_2022<(tfel::material::ModellingHypothesis::Hypothesis)6u, Type, false>::integrate(tfel::material::WoodOrthotropicElasticity_2022<(tfel::material::ModellingHypothesis::Hypothesis)6u, Type, false>::SMFlag, tfel::material::WoodOrthotropicElasticity_2022<(tfel::material::ModellingHypothesis::Hypothesis)6u, Type, false>::SMType)’: fort.5:63:28: error: ‘dw’ was not declared in this scope compilation terminated due to -Wfatal-errors. Makefile.mfront:37: recipe for target 'asterWoodOrthotropicElasticity_2022.o' failed make: *** [asterWoodOrthotropicElasticity_2022.o] Error 1 callMake: can't build target 'all' libraries building went wrong

thelfer commented 2 years ago

@lore-rip Did you turn w into a material property as in the initial implementation ? If so, this would explain the error as material properties are assumed to be given at the end of the time step.

lore-rip commented 2 years ago

@thelfer, you were right I have some other problems inside the code aster creation of the mfront library

CREA_LIB_MFRONT(DEBUG='NON', UNITE_LIBRAIRIE=7, UNITE_MFRONT=5)

Treating target : all The following library has been built :

Mémoire (Mo) : 987.58 / 979.47 / 47.79 / 24.37 (VmPeak / VmSize / Optimum / Minimum)

Fin commande #0003 user+syst: 2.05s (syst: 0.23s, elaps: 2.28s)

----------------------------------------------------------------------------------------------

Traceback (most recent call last): File "./point1.comm.changed.py", line 20, in UNITE_MFRONT=5) File "/opt/salome_meca/Salome-V2021-s9/tools/Code_aster_stable-1540/lib/aster/codeaster/Supervis/ExecuteCommand.py", line 172, in run return cmd.run(kwargs) File "/opt/salome_meca/Salome-V2021-s9/tools/Code_aster_stable-1540/lib/aster/codeaster/Supervis/ExecuteCommand.py", line 223, in run self.exec_(keywords) File "/opt/salome_meca/Salome-V2021-s9/tools/Code_aster_stable-1540/lib/aster/codeaster/Supervis/ExecuteCommand.py", line 703, in exec output = self._op(self, keywords) File "/opt/salome_meca/Salome-V2021-s9/tools/Code_aster_stable-1540/lib/aster/code_aster/MacroCommands/crea_lib_mfront_ops.py", line 50, in crea_lib_mfront_ops UTMESS('F', 'MFRONT_4', valk="libAsterBehaviour.so") File "/opt/salome_meca/Salome-V2021-s9/tools/Code_aster_stable-1540/lib/aster/code_aster/Messages/Utmess.py", line 653, in UTMESS exception=True, print_as=print_as, files=files, cc=cc) File "/opt/salome_meca/Salome-V2021-s9/tools/Code_aster_stable-1540/lib/aster/code_aster/Messages/Utmess.py", line 127, in call self.print_message(*args, **kwargs) File "/opt/salome_meca/Salome-V2021-s9/tools/Code_aster_stable-1540/lib/aster/code_aster/Messages/Utmess.py", line 188, in print_message raise libaster.AsterError(id0, valk, vali, valr) libaster.AsterError:

lore-rip commented 2 years ago

Anyways in my code aster full implementation the moisture field and the angle field are calculated a priori as NEUT fields. This assures the possibility also of using syrthes if needed. Attached the file I'm using roto_test1.mfront.zip

thelfer commented 2 years ago

@lore-rip you shall probably use

swell(0) = 0.00064 (NEUT2 + dNEUT2); swell(1) = 0.00115 (NEUT2 + dNEUT2); swell(2) = 0.00013;

BTW, what is the meaning of having a constant axial swelling ?

LoreRiparbelli commented 2 years ago

@thelfer it is clearly an error, sorry! BTW I'm coming with a very stupid question

you shall probably use

it make perfectly sense! Thanks

thelfer commented 2 years ago

Anyways in my code aster full implementation the moisture field and the angle field are calculated a priori as NEUT fields. This assures the possibility also of using syrthes if needed.

Are NEUT* special names ?

BTW, relative humidity is associated with the special name SECH (see https://thelfer.github.io/MFrontGallery/web/Burger_EDF_CIWAP_2021.html#variables-declarations). I imagine that MoistureContent and RelativeHumidity are related (egal ?) but did not dig further.

LoreRiparbelli commented 2 years ago

Yes Thomas. NEUT1 and NEUT2 are numerical fields calculated in code_aster. I don't use SECH for many reasons, the most important one is that I want to manipulate those fields, as you said to go from Relative Humidity to Moisture Content. The second one is the use of Syrthes. More generally to maintain a certain freedom from certain code_aster internal procedures

thelfer commented 2 years ago

@lore-rip, Indeed, code_aster is not very flexible when it comes to external state variables. Hopefully, we are discussing with the code_aster development team to improve that.

lore-rip commented 2 years ago

@thelfer Hi, this entry at the beginning of the script make the aster command CREA_LIB_MFRONT crash. @Material Wood; Maybe it is a bug.... Only to inform you why I was not able to run your script

lore-rip commented 2 years ago

Now your script works greatly, and with the following modifications it is possible to use it within SIMU_POINT_MAT in code_aster WoodOrthotropicElasticity_2022_a.mfront.zip

lore-rip commented 2 years ago

Summarizing we have 2 behaviors to describe cylindrical elasticity with swelling for wood:

1.WoodOrthotropicElasticity_2022_a.mfront that is the basic implementation of orthotropic behavior in cylindrical coordinates

2.WoodOrthotropicConicElasticity_2022_a.mfront that takes into account also of the possible conical growth of the rings in the vertical direction (here attached) WoodOrthotropicConicElasticity_2022_a.mfront.zip

lore-rip commented 2 years ago

@thelfer Hi, please find attached some regression test -Test_01a, Test_02a, Test_03a, are done WoodOrthotropicElasticity_2022_a.mfront for three angles 0°, 45°, 90° -Test_01b, Test_02b, Test_03b, are done WoodOrthotropicConicElasticity_2022_a.mfront for three angles 0°, 45°, 90°, and conical angle of 2° -Test_04b, Test_05b, Test_06b, are done WoodOrthotropicConicElasticity_2022_a.mfront for three angles 0°, 45°, 90°, and conical angle of 1°

Here the coefficient and the external variables

Test_01a MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 113., 1.0, 0.0, 0., 0.035))) Moisture content NEUT1 =1 Angle NEUT2 =0

Test_02a MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 113., 1.0, 0.0, 0., 0.035))) Moisture content NEUT1 =1 Angle NEUT2 =pi/4

Test_03a MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 113., 1.0, 0.0, 0., 0.035))) Moisture content NEUT1 =1 Angle NEUT2 =pi/2

Test_01b MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 113., 0.00064, 0.00013, 0.00115, 0.035))) Moisture content NEUT1 =1 Angle NEUT2 =0

Test_02b MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 113., 0.00064, 0.00013, 0.00115, 0.035))) Moisture content NEUT1 =1 Angle NEUT2 =pi/4

Test_03b MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 113., 0.00064, 0.00013, 0.00115, 0.035))) Moisture content NEUT1 =1 Angle NEUT2 =pi/2

Test_04b MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 113., 0.00064, 0.00013, 0.00115, 0.017))) Moisture content NEUT1 =1 Angle NEUT2 =0

Test_05b MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 113., 0.00064, 0.00013, 0.00115, 0.017))) Moisture content NEUT1 =1 Angle NEUT2 =pi/4

Test_06b MFRONT=_F(LISTE_COEF=(635., 10285., 369., 0.029, 0.42, 0.165, 786., 838., 113., 0.00064, 0.00013, 0.00115, 0.017))) Moisture content NEUT1 =1 Angle NEUT2 =pi/2 tests_mfront.zip