Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Wood material repo #2

Open
lore-rip opened this issue Dec 9, 2021 · 41 comments
Open

Wood material repo #2

lore-rip opened this issue Dec 9, 2021 · 41 comments

Comments

@lore-rip
Copy link

lore-rip commented Dec 9, 2021

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
Copy link
Owner

thelfer commented Dec 10, 2021

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
Copy link
Author

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
Copy link
Owner

thelfer commented Dec 16, 2021

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

@thelfer
Copy link
Owner

thelfer commented Dec 16, 2021

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:

  • The name of the behaviour dedicated to poplar. 2021 is an identifier which shall probably be more specific and reasonnably describe this version of the behaviour. If you did the idenfification, something like FlorenceUniversity2021 would do.
  • Orthotropic axes conventions which may require that we change the order of the parameters. See this paragraph to have some insights about this question http://tfel.sourceforge.net/faq.html#orthotropic-axes-convention. I would recommend to follow the Pipe orthotropic convention.

@lore-rip
Copy link
Author

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
Copy link
Author

lore-rip commented Dec 17, 2021

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
Copy link
Owner

thelfer commented Dec 17, 2021

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
Copy link
Author

lore-rip commented Dec 17, 2021

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

@thelfer
Copy link
Owner

thelfer commented Dec 17, 2021

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
Copy link
Owner

thelfer commented Dec 18, 2021

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
Copy link
Author

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

@thelfer
Copy link
Owner

thelfer commented Dec 18, 2021

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

@lore-rip
Copy link
Author

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

@thelfer
Copy link
Owner

thelfer commented Dec 20, 2021

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
Copy link
Author

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
Copy link
Owner

thelfer commented Mar 31, 2022

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

@thelfer
Copy link
Owner

thelfer commented Apr 1, 2022

@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
Copy link

LoreRiparbelli commented Apr 1, 2022

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
Copy link
Author

lore-rip commented Apr 1, 2022

@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
Copy link
Owner

thelfer commented Apr 1, 2022

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

@LoreRiparbelli
Copy link

@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
Copy link
Owner

thelfer commented Apr 1, 2022

@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
Copy link

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
Copy link
Owner

thelfer commented Apr 1, 2022

You're right. I got your mail.

@LoreRiparbelli
Copy link

Thank you very much!! :)

thelfer added a commit that referenced this issue Apr 1, 2022
@thelfer
Copy link
Owner

thelfer commented Apr 1, 2022

@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
Copy link
Author

lore-rip commented Apr 2, 2022

@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
Copy link
Owner

thelfer commented Apr 2, 2022

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

@lore-rip
Copy link
Author

lore-rip commented Apr 2, 2022

@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
Copy link
Owner

thelfer commented Apr 2, 2022

@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
Copy link
Author

lore-rip commented Apr 3, 2022

@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 :

  • libAsterWood.so : asterwoodorthotropicelasticity_2022

╔════════════════════════════════════════════════════════════════════════════════════════════════╗
║ <MFRONT_4> ║
║ ║
║ Le fichier de sortie de MFront libAsterBehaviour.so n'a pas pu être produit. ║
║ ║
║ Conseil : Vérifiez l'existence du fichier fourni au mot-clé UNITE_MFRONT. ║
╚════════════════════════════════════════════════════════════════════════════════════════════════╝

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

Fin commande #3 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/code_aster/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/code_aster/Supervis/ExecuteCommand.py", line 223, in run_
self.exec_(keywords)
File "/opt/salome_meca/Salome-V2021-s9/tools/Code_aster_stable-1540/lib/aster/code_aster/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
Copy link
Author

lore-rip commented Apr 3, 2022

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
Copy link
Owner

thelfer commented Apr 3, 2022

@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
Copy link

@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
Copy link
Owner

thelfer commented Apr 4, 2022

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
Copy link

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
Copy link
Owner

thelfer commented Apr 4, 2022

@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
Copy link
Author

lore-rip commented Apr 5, 2022

@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
Copy link
Author

lore-rip commented Apr 5, 2022

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
Copy link
Author

lore-rip commented Apr 5, 2022

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
Copy link
Author

lore-rip commented Apr 7, 2022

@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

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

3 participants