argonne-national-laboratory / neml Goto Github PK

In Blackbear we need to be able to specify parameters in the MOSOE input file (this is required for using the stochastic methods module to launch many simulations). We currently use a hack, where a {variable} string in the XML code is replaced through string operation in a preprocessing step.

The XML spec has the concept of "entities" (&VARIABLE;) which unfortunately is not supported by the rapidxml library, so we need to stick with find/replace for now.

However, I'd like to perform that find/replace on individual attributes in the DOM tree. This would permit better error checking.

Right now, in Blackbear we cannot do proper error checking for missing or unused substitutions, because we cannot limit the find/replace to only the model subblock that is actually used.

This is weird: I can configure cmake and compile just fine without the python binding. But if I set WRAP_PYTHON=ON I get the following error:

CMake Error in src/CMakeLists.txt:
  The compiler feature "cxx_inheriting_constructors" is not known to CXX
  compiler

  "GNU"

  version 9.3.1.

However, even after the error, cmake is able to write build files, and it seems that neml compiles fine.

We've found an odd corner case where the linesearch system fails to find a good step, even though if you take a random small step in the direction of the NR increment everything works fine on subsequent iterations. I can't reproduce it with a simple case, but once we can we should figure out what's going on.

Hello,
Tried building neml on Win10 machine with installed LAPACK and BLAS through OpenBlas, but cmake invoke gives me an error that it is unable to find BLAS. Has anybody experienced similar issue ?:

D:\GITlib\neml>cmake -DCMAKE_CXX_COMPILER=icl.exe -DBLA_VENDOR=OpenBLAS -DBoost_INCLUDE_DIR="D:\GITlib\boost" -DCMAKE_BUILD_TYPE=Release -DUSE_OPENMP=False -G Ninja .
CMake Error at D:/GITlib/cmake/share/cmake-3.20/Modules/FindPackageHandleStandardArgs.cmake:230 (message):
  Could NOT find BLAS (missing: BLAS_LIBRARIES)
Call Stack (most recent call first):
  D:/GITlib/cmake/share/cmake-3.20/Modules/FindPackageHandleStandardArgs.cmake:594 (_FPHSA_FAILURE_MESSAGE)
  D:/GITlib/cmake/share/cmake-3.20/Modules/FindBLAS.cmake:1045 (find_package_handle_standard_args)
  CMakeLists.txt:42 (FIND_PACKAGE)


-- Configuring incomplete, errors occurred!
See also "D:/GITlib/neml/CMakeFiles/CMakeOutput.log".
See also "D:/GITlib/neml/CMakeFiles/CMakeError.log".

A second point - neml looks like a great suite of models, but it appears the community can only use \util for Abaqus. Is there an equivalent for ANSYS? Thx

Add an auxkernel computing traction on a boundary.
The user can select to output either a traction component, or the norm of the shear traction.

Huge project: add single crystal plasticity to NEML.

The following file has a typo on line 9.
<plastic type = type="GeneralIntegrator">
GR91CPfit.txt

The parser errors out but it does not provide a line number, which it should.
Andrea

The recent moving of the headers from neml/src to neml/include broke our build of blackbear against the current neml. That was easy to fix by using neml/include as our include directory. However, another change was made within neml that makes it so we still can't build. There used to be a #include "models.h" at the top of neml_interface.h, without which we can't build against NEML because NEMLModel is undefined. I think we just need to add that back in.

Add LLNL's SNLS solver to compare to our built-in Newton and SNL's NOX.

NEML was setup to provide several stress update methods for each constitutive model. Up to now we have just been implementing small strain models and provided only a small strain material update. We should support, at least, an incremental large strain formulation using an objective rate to take the small strain stress rate and consistently apply it to large strain kinematics.

For continuation and recovery it would be nice if the internal state of NEML models could be serialized and unserialized.

NEML has a pattern of passing up ier return status values up the call tree. This pattern is broken in solvers.cxx where teh return value of system->RJ is ignores in multiple places.

I'd like to return a new error status (something like domain error or whatever you'd like to call it) if results are outside of a reasonable range of values and would result in fatal errors further down the line.

This method is required have a user friendly approach fro retryving history variables by names in any linked software.
For instance in Deer we store history variables required by neml, but we can't acces them by nam,e because deer doesn't know which location on the history vector is assigned to a particular history variable.
This make our life much easier for postprocessing

In objects.h is a bit of code that looks like it would run an infinite loop:

  /// Create and cast an object to a type as a unique_ptr
  template<typename T>
  std::unique_ptr<T> create_unique(ParameterSet & params)
  {
    auto res = std::unique_ptr<T>(dynamic_cast<T*>(create_unique(params).release()));
    if (res == nullptr) {
      throw WrongTypeError();
    }
    else {
      return res;
    }
  }

The function calls itself. It seems like the first line should be

    auto res = std::unique_ptr<T>(dynamic_cast<T*>(create(params).release()));

To execute the shared pointer version. I'm a bit confused as to why this code is never hit though (I probably took a wrong turn tracing the call chain) :-)

Will it be a good idea if we pass the history not as a history variable but a double const * h_np1 ( similar to set_active_orientation). That way we will, most likely, have a consistent way passing the history between DEER and NEML. For example the strength function being:

double SingleCrystalModel::strength(double * const hist, Lattice & L,
double T) const
{
History h = gather_history_(hist);
return kinematics_->strength(h, L, T);
}

Installation on Linux was super straight forward, but on mac we ran into a some small bumps. I have prepared a small MD document to document the prerequisites needed (for mac homebrew and ubuntu apt)

It would be neat if NEML could provide some feedback to the embedding code about the suitability of the current timestep.

My problems reported in the other issue basically vanished when I manually cut the timestep at a critical point in the simulation. I was ramping up the temperature and let MOOSE increase the timestep as the solve converged very well, but at a certain temperature the simulation went from converging very well to failing during NEML stress updates. There is not indicator that MOOSE can use to prevent a runaway increase in timestep to avoid the NEML failures in the first place. I'm still trying to figure out why MOOSE doesn't even cut the timestep reliable when the NEML failure occurs, but that's an orthogonal problem.

Is there a way to query the number of internal NEML iterations maybe?

Some pointless change that someone wants.

We were trying to add a couple of simple regression tests that compare against native MOOSE models for linear strain hardening plasticity and power law creep in this BlackBear PR:
idaholab/blackbear#65
The power law creep model encounters an out-of-bounds access error on line 1088 of models.cxx when run with the debug options that we use because it's trying to access index 0 of a history variable array that happens to be of size 0 in our case because we're using a perfect plasticity model.

I think the right fix here is to check for that case and pass in a null pointer to the function that needs the history data. I'll submit a PR with a proposed fix.

Const correctness is violate with new libxml++3 versions. The get_children method on a const xmlpp::Node returns a list of const xmlpp::Node * (so that whole method needs to return a const xmlpp::Node *.

I have a PR ready...

It would be nice to have a built-in newton solver that can take advantage of having "the action of the inverse jacobian" rather than the jacobian matrix itself.

There are two problem with how we implement damage:

1. NEMLScalarDamageModel_sd is way too tightly coupled. This class should take as input a base model and a scalar damage model implemented as a separate class. The creep models are a good example of what this should look like.

2. I can't for the life of me remember why I implemented these in the discretized form where you're expected to directly mess with the stress update, instead of in the rate form like all of the rest of the code. So you should provide a damage evolution rate and the superclass should take care of the integration. I have some vague recollection of rejecting this approach for a good reason, but I sure can't think of it now.

These BLAS/Lapack forward declarations conflict with declaration we pull in through PETSc:

See compile error here: https://civet.inl.gov/job/2140682/

Compiling C++ (in opt mode) /data/civet0/build/blackbear/moose/framework/build/unity_src/distributions_Unity.C...
In file included from /data/civet0/build/blackbear/contrib/neml/include/math/rotations.h:5,
                 from /data/civet0/build/blackbear/contrib/neml/include/elasticity.h:7,
                 from /data/civet0/build/blackbear/contrib/neml/include/models.h:5,
                 from /data/civet0/build/blackbear/contrib/neml/include/parse.h:5,
                 from /data/civet0/build/blackbear/moose/framework/src/utils/InputParameters.C:24,
                 from /data/civet0/build/blackbear/moose/framework/build/unity_src/utils_Unity.C:28:
/data/civet0/build/blackbear/contrib/neml/include/math/nemlmath.h:19:8: error: conflicting declaration of C function 'void dgetrf_(const int&, const int&, double*, const int&, int*, int&)'
   19 |   void dgetrf_(const int & m, const int & n, double* A, const int & lda, int* ipiv, int & info);
      |        ^~~~~~~
In file included from /opt/petsc/include/petscblaslapack.h:72,
                 from /data/civet0/build/blackbear/moose/framework/src/utils/ColumnMajorMatrix.C:19,
                 from /data/civet0/build/blackbear/moose/framework/build/unity_src/utils_Unity.C:9:
/opt/petsc/include/petscblaslapack_mangle.h:53:47: note: previous declaration 'void dgetrf_(const PetscBLASInt*, const PetscBLASInt*, PetscReal*, const PetscBLASInt*, PetscBLASInt*, PetscBLASInt*)'

Ping @bwspenc

I want to have a generic work-based damage model for creep-fatigue in NEML.

Looks like using an enum as an unordered_map index is not fully portable. This can be mitigated by passing an int hasher (std::hash<int>) as a third template argument to the unordered map type.

A scaling prefactor is useful when converting material parameters from macroscopic to microscopic

@reverendbedford @hugary1995
I am using NEML within MOOSE, and I am trying to output the damage index for visualization.
NEMLScalarDamagedModel_sd::get_damage is returning the wrong answer. Its seems to be returning the the value at the wrong index.
I'll try to explain as much as I can, using debugging outputs that I have. The problem I am running should have damage everywhere within the first time step.
The history as a quadrature point is below. Given that I know there is damage and that the damage should be between 0 and 1, I assume that the last value is the correct value.

_hist[_qp]
size=8
[0]: 0
[1]: 0
[2]: 0
[3]: 0
[4]: 0
[5]: 0
[6]: 10.025091709671862
[7]: 0.47206766361193314

Looking in NEMLScalarDamagedModel_sd::get_damage(const double *const h_np1)

History h = gather_history_(h_np1);
return h.get<double>(prefix("damage"));

we can see that the actual return value will be the 6th index which would give an answer of 10.

h.loc_
size=3
[0]: {first:"small_stress", second:0}
[1]: {first:"damage", second:6}
[2]: {first:"alpha", second:7}

The XML file is below, I only added the <alpha> 0.99 </alpha> parameter for a hunch that alpha is never outputted into history but messes up the indexing, but I don't know if its supposed to actually be placed into history anyway. Let me if you need more information.

<materials>
  <test_powerdamage type="NEMLScalarDamagedModel_sd">
       <dkill> 0.25</dkill>
       <alpha> 0.99 </alpha>
    <elastic type="IsotropicLinearElasticModel">
      <m1>92000.0</m1>
      <m1_type>youngs</m1_type>
      <m2>0.3</m2>
      <m2_type>poissons</m2_type>
    </elastic>

    <base type="SmallStrainRateIndependentPlasticity">
      <elastic type="IsotropicLinearElasticModel">
        <m1>92000.0</m1>
        <m1_type>youngs</m1_type>
        <m2>0.3</m2>
        <m2_type>poissons</m2_type>
      </elastic>
      <flow type="RateIndependentAssociativeFlow">
        <surface type="IsoJ2"/>
        <hardening type="LinearIsotropicHardeningRule">
          <s0>180.0</s0>
          <K>1000.0</K>
        </hardening>
      </flow>
    </base>

    <damage type="PowerLawDamage">
      <elastic type="IsotropicLinearElasticModel">
        <m1>92000.0</m1>
        <m1_type>youngs</m1_type>
        <m2>0.3</m2>
        <m2_type>poissons</m2_type>
      </elastic>
      <A>2.0e-5</A>
      <a>2.2</a>
    </damage>

  </test_powerdamage>

</materials>

Right now the error system is a mess. It's error codes in most of the code which are checked and turned into exceptions for the python bindings. But the object system and the parser throws exceptions, which are caught and turned into codes for the C and fortran bindings.

I would like to move over to exceptions in the whole codebase (catching them before going into the C and fortran interfaces). This will involve:

1. Lots of busy work designing appropriate exceptions and changing the function signature in nearly every method and function in NEML.
2. Going into adaptive integration loops and changing error code checks into catches.
3. Registering the errors with pybind11 and making sure they translate correctly in the python bindings.

Do all the hooks work? Make a pointless change!

I'm having some difficulties debugging NEML runs where I run into floating point exceptions. The problem here is that the actual FPE comes fairly late in the game (eigenvalue decomposition), after previous calculations returned Inf several steps before in a separate branch of the call tree (that is not in the backtrace).

In my particular case the first accident is occurring in the hardening rule:

int VoceIsotropicHardeningRule::q(const double * const alpha,
                                    double T, double * const qv) const
{
  qv[0] = -s0_->value(T) - R_->value(T) * (1.0 - exp(-d_->value(T) * alpha[0]));
 ...

Where qv[0] ends up being Inf because alpha[0] < -1000.0 and d_->value(T) ~ 1.

I suggest adding a bunch of assertions... and maybe a few more comments in the code to explain teh physical meaning of each calculated quantity...

Right now we pass in a huge list of parameters to update the material model state. This is both ugly and makes it very hard to extend to add new types of state descriptions (e.g. radiation dose).

We have the same problem when it comes to history variables. Right now we just pass around vectors of doubles. We need a way to label variables with names and accommodate non scalar state variables.

I think we could deal with both of these issues with the same data structure. Requirements:

1. Store scalars, vectors, symmetric 2nd order tensors (Mandel notation), and general 2nd order tensors along with an instance of this same data structure (to store history variables).
2. Be able to report which type each of the above has.
3. Associate a name with each variable and be able to retrieve variables by name.
4. Be able to marshal the struct into a flat vector and reconstruct the structure from the marshaled vector.

Objects would all need to be updated to:

1. Report what sort of state they need (both history and general state).
2. Take the state structs instead of flat vectors as input.

Interfaces would all need to be updated to:

1. Provide the struct instead of enumerated parameters (or a flat vector for history).

We are having trouble providing the libxml++ library on all platforms of interest.

I propose switching over to use a different XML parser. RapidXML (http://rapidxml.sourceforge.net/) looks like a good option because it's header-only and has a permissive license, so we can just bundle it with NEML in a contrib directory and #include the header files to use it.

We want a block evaluation routine that takes tensors as input, coverts to Mandel efficient, does the stress update for the whole block, and efficiently reconverts to tensors.

Time to move to python3. 2.7 will lose support in a few months.

People want to compile and run NEML on Windows, particularly as a UMAT.

For unknown reasons some build platforms don't get the cmath header in larsonmiller. Can be fixed by explicitly adding an "include nemlmath".

In certain models the value for yv in general_flow.cxx:81 can get very large. This will cause a LAPACK call in line 103 to throw a floating point exception.

In MOOSE this hardware exception will be caught by libMesh and result in the termination of the simulation. What we actually want is soft abort the current step and retry with a smaller time step. If I insert a line 83

  if (yv > 1e10) return LINALG_FAILURE;

This works beautifully. I picked 1e10 randomly (1e50 still led to the FPE), which is not quite satisfactory. If you have more physical insight you may have a better rationale for a threshold value.

I'm trying to compile NEML as a library as part of the MOOSE build process, rather than using cmake, and running into errors like this:

neml/src/surfaces.cxx:473:10: error: call to 'make_unique' is ambiguous

I find myself wanting a "generic" creep rate model that returns the log of the creep rate as a generic Interpolate function of the log of the stress.

So the model will be: log(creep rate) = f(log(stress)) where f is one of our interpolate functions.

Add the MeanEffectiveStress defined as tr(sigma)/3 and fix HudllestonEffectiveStress to work when ss=0

Compiling BlackBear (which includes NEML) in debug mode shows this warning that indicates a potentially serious issue:

/Users/reverendbedford/projects/blackbear/contrib/neml/src/objects.h:121:12: warning: returning address of local temporary object [-Wreturn-stack-address]
    return ss.str().c_str();
           ^~~~~~~~

(User name swapped 😄)

NEML (at least the parts I use, not the wrappers) depends on the Boost library because of its use of variant. Although we (MOOSE) support pulling in the full Boost library, it's a non-standard way of building the code, and has been causing an unreal amount of headaches with our testing system because of unintended side effects.

To give an example, we were getting inundated with compiler warning messages coming from another boost file we were including (unrelated to variant), which actually made our builds fail because the system couldn't handle the sheer volume of output.

At a minimum, NEML's dependency on boost is a barrier to its more widespread use because using it requires recompiling libmesh with boost, which is nonstandard and enough of a pain to do that someone has to be serious about using those models to go to that trouble.

Anyway, suffice it to say that I think that the small benefit obtained from using variant doesn't justify the pain of building with boost. I see that variant is supported in C++17, which we are talking about moving to, so that could solve that problem, but we're a little ways off from being there. I would propose that NEML use a design more like what we use in libmesh/MOOSE, where we have templated classes that store parameters of various types. It looks to me like the effects of that change would be pretty localized to a small part of the NEML codebase.

Since I'm asking for this, I'm happy to implement it. If you guys absolutely hate what I come up with, maybe you can move back to using variant in a while when everyone is using C++17.

For integration with MOOSE it would be nice to be able to create NEML classes templated on either regular floating point numbers as well as dual numbers (to compute derivatives using forward mode automatic differentiation).

I guess this would make NEML a largely header only library which might no be ideal for other users though :-/

Just throwing this out.

Our test suite only checks things through the bindings and never checks the standalone library. We should test both.