Multi Scalar Field Element Development

[Roozbeh Rezakhani]

Hello Everyone,

I have successfully developed 1D element for the diffusion of a scalar field (temperature) using "ChNodeFEAxyzP" in /src/chrono/fea.

Next, I generalized the 1D routine for the diffusion of 3 scalar fields (temperature, humidity, chemical content). After compilation, I get the following error during the first solve step:

FEA_hydrothermal_2D_square: /usr/include/eigen3/Eigen/src/Core/Block.h:146: Eigen::Block<XprType, BlockRows, BlockCols, InnerPanel>::Block(XprType&, Eigen::Index, Eigen::Index, Eigen::Index, Eigen::Index) [with XprType = const Eigen::Matrix<double, -1, 1>; int BlockRows = -1; int BlockCols = 1; bool InnerPanel = false; Eigen::Index = long int]: Assertion `startRow >= 0 && blockRows >= 0 && startRow <= xpr.rows() - blockRows && startCol >= 0 && blockCols >= 0 && startCol <= xpr.cols() - blockCols' failed.
Aborted (core dumped)

I would be thankful for your input.

Thanks,

Roozbeh

[Luning Bakke]

Hello Roozbeh,

While my background is not solving thermal PDEs, I can provide some ideas. I also need some clarifications. 

• Are the three scalar fields (temperature, humidity, and chemical content) coupled or decoupled? Implementing a solver in Chrono with all three fields coupled is a lot more complicated, and I suggest starting with a small problem with 2 fields coupled. 
• When you mention that you have "generalized the 1D routine for the diffusion of 3 scalar fields", does it mean, 
  • you customized your own node class based on the ChNodeFEAxyzP.h to handle 3 scalar fields? That's doable, but you have to properly initialize, populate and assemble all the fields. 
  • Or you created three separate meshes, each for one field variable? This should work for a decoupled case; make sure you adjust the thermal conductivity and specific heat values to match the governing equations for humidity and chemical content. You can still use the same API call, since under the hood, the same type of equations is solved. 
• Regarding the the error you saw from Eigen, this usually happens when some matrix/vector access is out of bounds. Maybe the indices are incorrect, or the fields are not initialized properly. I suggest that you step through the code with a debugger, check the dimensions, and also verify the matrix/vector size after initialization. 

Hope my suggestions help.

Thank you,

Luning