Optcode

doc/elementlibmanual/elementlibmanual.tex

@@ -437,7 +437,7 @@ \subsubsection{Lattice3dBoundary element}
 \elementDescription{Reference}{\cite{AthWheGra18}}
 \end{elementsummary}
 
-\subsubsection{Latticelink3d}
+\subsubsection{Latticelink3d element}
 This element represents a two-node 3d link element connecting 3d beam and 3d lattice elements. Each node has six degrees of freedom. 
 The input parameters for this element are shown in Table~\ref{latticelink3dsummary}. 
 
@@ -451,14 +451,41 @@ \subsubsection{Latticelink3d}
 \elementDescription{Reference}{\cite{GraAnt19}}
 \end{elementsummary}
 
-\subsubsection{Latticelink3dboundary}
+\subsubsection{Latticelink3dboundary element}
 Represents three-node 3d boundary link element connecting 3d beam and 3d lattice elements. The first two nodes have the same meaning as for latticelink3d. The third node is used to control the loading of the periodic cell using three normal (xx, yy and zz) and three shear strain (yz, xz, xy) components.
 The specific input parameters for this element in addition of those for latticelink3d are shown in Table~\ref{latticelink3dboundarysummary}. 
 
 \begin{elementsummary}{latticelink3dboundary}{3d lattice link boundary element}{{ \field{location}{in} }{rn}}{latticelink3dboundary element summary}{latticelink3dboundarysummary}
 \elementParam{\param{location}: }
 \elementDescription{Unknowns}{Six dofs ($u$-displacement, $v$-displacement, $w$-displacement, $u$-rotation, $v$-rotation and $w$-rotation) are required in each node.}
 \elementDescription{Reference}{\cite{GraAnt19}}
+\end{elementsummary}
+
+ \subsubsection{Latticeframe3d element}
+ Latticeframe3d represents a two-node linear 3d frame element. Each node has six degrees of freedom, namely three translations and three rotations. The element is based on the link between rigid body spring and Timoshenko frame models. The main idea is to model the elastic and inelastic response of a connection of two nodes by a set of springs located at the contact facet of two rigid bodies. Displacement jumps are computed at the contact facet, which are smeared out over the element length in the form of strains. The input parameters for this element are shown in Table ~\ref{latticeframe3dsummary}. For the computation of the nodal forces, the initial configuration is used.
+%%  \begin{figure}[htb]
+%%  \centering
+%%  \includegraphics[width=0.3\textwidth]{./latticeframe.png}
+%%  \includegraphics[width=0.3\textwidth]{./latticeframe1.png}
+%%  \caption{latticeframe3d element. Node numbering, DOF numbering, cross-section vertices and local coordinate system at integration point.}
+%%  \label{latticeframe3dfig}
+%%  \end{figure}
+
+ \begin{elementsummary}{latticeframe3d}{3d lattice frame element}{{ \field{zaxis}{ra} }{rn}}{latticeframe3d element summary}{latticeframe3dsummary}
+\elementParam{\param{zaxis}: z-axis of local coordinate system }
+\elementDescription{CS properties}{Area, inertia moment along y and z axis (\param{iy} and \param{iz} parameters), torsion inertia moment (\param{ik} parameter) and either cross section area shear correction factor (\param{beamshearcoeff} parameter) or equivalent shear areas (\param{shearareay} and \param{shearareaz} parameters) are required. These cross section properties are assumed to be defined in local coordinate system of element.}
+\elementDescription{Unknowns}{Six dofs ($u$-displacement, $v$-displacement, $w$-displacement, $u$-rotation, $v$-rotation and $w$-rotation) are required in each node.}
+\elementDescription{Reference}{\cite{Toi91, Toi93}}
+\end{elementsummary}
+
+ \subsubsection{Latticeframe3dnl element}
+ Latticeframe3dnl represents a two-node nonlinear 3d frame element. Each node has six degrees of freedom, namely three translations and three rotations. The element is based on the link between rigid body spring and Timoshenko frame models. The main idea is to model the elastic and inelastic response of a connection of two nodes by a set of springs located at the contact facet of two rigid bodies. Displacement jumps are computed at the contact facet considering the incrementally deformed geometry, which are smeared out over the element length in the form of strains. The input parameters for this element are shown in Table ~\ref{latticeframe3dnlsummary}. For the computation of the nodal forces, the incrementally deformed configuration is used. At the start of each step the local co-ordinate system is updated with the displacements of the previous step.
+
+ \begin{elementsummary}{latticeframe3dnl}{3d lattice frame element}{{ \field{zaxis}{ra} }{rn}}{latticeframe3d element summary}{latticeframe3dnlsummary}
+\elementParam{\param{zaxis}: z-axis of local coordinate system }
+\elementDescription{CS properties}{Area, inertia moment along y and z axis (\param{iy} and \param{iz} parameters), torsion inertia moment (\param{ik} parameter) and either cross section area shear correction factor (\param{beamshearcoeff} parameter) or equivalent shear areas (\param{shearareay} and \param{shearareaz} parameters) are required. These cross section properties are assumed to be defined in local coordinate system of element.}
+\elementDescription{Unknowns}{Six dofs ($u$-displacement, $v$-displacement, $w$-displacement, $u$-rotation, $v$-rotation and $w$-rotation) are required in each node.}
+\elementDescription{Reference}{\cite{Toi91, Toi93, AbdGra23}}
 \end{elementsummary}
 
 %-----------------------------------------------------------------------------------------------
@@ -2240,17 +2267,15 @@ \subsubsection{Tet1\_3D\_SUPG element}
 %\printbibliography
 
 \begin{thebibliography}{99}
-
-
 \bibitem{RobertCook1989} R. D. Cook and D. S. Malkus and M. E. Plesha, ``Concepts and Applications of Finite Element Analysis'', Third Edition, isbn: 0-471-84788-7, 1989.
 \bibitem{BittnarSejnoha1996} Z. Bittnar and J. Sejnoha, ``Numerical Methods in Structural Mechanics'',Thomas Telford,isbn:978-0784401705, 1996.
 \bibitem{RagnarLarsson2011} R. Larsson and J. Mediavilla and M. Fagerström, ``Dynamic fracture modeling in shell structures based on XFEM'', International Journal for Numerical Methods in Engineering, vol. 86, no. 4-5, 499--527, 2011.
 \bibitem{GraJir10} P. Grassl and M. Jir\'{a}sek, ``Meso-scale approach to modelling the fracture process zone of concrete subjected to uniaxial tension'', International Journal of Solids and Structures, vol. 47, iss. 7-8, pp. 957-968, 2010..
 \bibitem{GraBol16} P. Grassl, J. Bolander, ``Three-Dimensional Network Model for Coupling of Fracture and Mass Transport in Quasi-Brittle Geomaterials'', Materials, 9, 782, 2016
 \bibitem{AthWheGra18} I. Athanasiadis, S. Wheeler and P. Grassl. ``Hydro-mechanical network modelling of particulate composites'', International Journal of Solids and Structures, vol. 130-131, pp. 49-60, 2018.
-\bibitem{GraAnt19} P. Grassl and A. Antonelli. ``3D network modelling of fracture processes in fibre-reinforced geomaterials'', International Journal of Solids and Structures, vol. 156-157, Pages 234-242, 2019.
-
-  
+\bibitem{GraAnt19} P. Grassl and A. Antonelli. ``3D network modelling of fracture processes in fibre-reinforced geomaterials'', International Journal of Solids and Structures, vol. 156-157, pp. 234-242, 2019.
+\bibitem{Toi91} Y. Toi. ``Shifted integration technique in one‐dimensional plastic collapse analysis using linear and cubic finite elements'', International Journal for Numerical Methods in Engineering 31, no. 8, pp. 1537-1552, 1991.
+\bibitem{AbdGra23} G. Abdelrhim and P. Grassl. ``On a geometric nonlinear Timoshenko frame element with material nonlinearity'', Annual Conference of the UK Association for Computational Mechanics (UKACM), The University of Warwick, UK, 19-21 April 2023.
 \end{thebibliography}
 \end{document}
 

doc/matlibmanual/matlibmanual.tex

@@ -6079,6 +6079,74 @@ \subsubsection{Viscoelastic plasticity damage lattice model}
 \label{latticeplasticdamageviscoelastic_table}
 \end{table}
 
+\subsubsection{Plasticity lattice model for steel}
+
+This is a Plasticity material used together with lattice elements.It has bean developed to describe the failure process  failure of steel frame member.
+The stress-strain law has the form
+\begin{equation}
+\boldsymbol{\sigma} = \mathbf{D}_{\rm e} \left(\boldsymbol{\varepsilon}-\boldsymbol{\varepsilon}_{\rm p} \right)
+\end{equation}
+where $\boldsymbol{\sigma}$ is a vector of tractions and rotational components, $\mathbf{D}_{\rm e}$ is the elastic stiffness matrix, $\boldsymbol{\varepsilon}$ is a vector of strains obtained from displacement jumps smeared over the element length and rotational components and $\boldsymbol{\varepsilon}_{\rm p}$ are the plastic strains. The model parameters are summarised in Table \ref{latticeframesteelplastic_table}.
+\begin{table}[!htb]
+ \centering
+\begin{tabular}{|l|l|}
+\hline
+Description & Plasticity lattice model for steel \\
+\hline
+Record Format & \descitem{latticeframesteelplastic} $\elemparam{}_{(in)}$
+$\elemparam{d}_{(rn)}$ $\optelemparam{talpha}_{(rn)}$ $\elemparam{e}_{(rn)}$ $\optelemparam{n}_{(rn)}$ $\optelemparam{nx0}_{(rn)}$ $\optelemparam{mx0}_{(rn)}$\\ &$\optelemparam{my0}_{(rn)}$ $\optelemparam{mz0}_{(rn)}$ $\elemparam{sub}_{(in)}$ $\elemparam{iter}_{(in)}$ $\elemparam{tol}_{(rn)}$\\
+Parameters &- \param{} material number\\
+&- $\param{d}$ material density\\
+&- $\optparam{talpha}$ Thermal expansion coefficient. Default is 0.\\
+&- $\optparam{e}$ Young's modulus of the lattice material.\\
+&- $\param{n}$ Poisson's ratio of the material that the beam element is made of.\\
+&- $\optparam{nx0}$ ultimate capacity under pure axialloads.\\
+&- $\optparam{mx0}$ ultimate capacity under pure moment about $x$.\\
+&- $\optparam{my0}$ ultimate capacity under pure moment about $y$.\\
+&- $\optparam{mz0}$ ultimate capacity under pure moment about $z$.\\
+&- $\optparam{sub}$ maximum number of subincrementations.\\
+&- $\optparam{iter}$ maximum number of newton iterations.\\
+&- $\optparam{tol}$ tolerance for newton method.\\
+Supported modes & 3dlattice\\
+
+\hline
+\end{tabular}
+
+\caption{Plasticity steel model for lattice elements -- summary.}
+\label{latticeframesteelplastic_table}
+%\end{center}                  
+  \end{table}
+The nonlinear response of the beam element is given directly as an relationship between internal force and strain. The basic equations include an additive decomposition of total strain into elastic part and plastic part
+
+\begin{equation}
+\boldsymbol{\varepsilon} = \boldsymbol{\varepsilon}_{\rm e}-\boldsymbol{\varepsilon}_{\rm p} 
+\end{equation}
+
+Trail stress
+\begin{equation}
+    \sigma_{n+1}^{tr}  =D^e (\varepsilon_{n+1}  - \varepsilon_n^p)
+\end{equation}
+
+The yield function
+\begin{equation}
+(\frac{N_x}{N_0})^2+ (\frac{M_x}{M_{x0}} )^2+ (\frac{M_y}{M_{y0}} )^2+ (\frac{M_z}{M_{z0}} )^2- 1 = 0
+\end{equation}
+
+where $M$, $M y$ , $N$ and $M$ , are the two components of bending moment, an axial force and
+a torsional moment, respectively. The subscript 0 indicates a fully plastic value under the
+condition that each component of resultant forces acts independently on a cross-section of
+a member.
+\\
+loading-unloading conditions
+
+loading-unloading conditions
+\begin{equation}
+f(\bar{\vsig},\kappa)\le 0 \qquad \dot{\lambda}\geq 0 \qquad \dot{\lambda}f(\bar{\vsig},\kappa)=0,
+\end{equation}
+evolution law for plastic strain
+\begin{equation}
+\dot{\veps}_{\rm{p}} = \dot{\lambda} \frac{\partial f}{\partial \bar{\vsig}},
+\end{equation}
 
 
 \subsection{Material models for steel relaxation}
@@ -6191,7 +6259,6 @@ \subsubsection{Model for relaxation of prestressing steel - SteelRelaxMat}
 \end{table}
 
 
-
 \subsection{User-defined material models using MFront}
 A user-defined material created with MFront can be used if MGIS is installed
 \newline

src/sm/CMakeLists.txt

@@ -97,6 +97,7 @@ set (sm_element
     Elements/LatticeElements/latticebeam3d.C
     Elements/LatticeElements/latticebeam3dboundary.C
     Elements/LatticeElements/latticeframe3d.C
+    Elements/LatticeElements/latticeframe3dnl.C
     Elements/tet21ghostsolid.C
     Elements/quad1platesubsoil.C 
     Elements/quad2platesubsoil.C
@@ -297,9 +298,9 @@ set (sm_material
     Materials/LatticeMaterials/latticeplasticitydamageviscoelastic.C
     Materials/LatticeMaterials/latticeframeelastic.C
     Materials/LatticeMaterials/latticeframesteelplastic.C
+    Materials/LatticeMaterials/latticeframeconcreteplastic.C
     Materials/1D/isoasymm1d.C
 )
-
 if (USE_MFRONT)
     list (APPEND sm_material Materials/mfrontusermaterial.C)
 endif ()
@@ -429,8 +430,9 @@ set (sm
     ${sm_obsolete}
     ${sm_new}
     ${sm_boundary_conditions}
-    ${sm_quasicontinuum} 
-    )
+    ${sm_quasicontinuum}
+        Elements/LatticeElements/latticebeam3dboundary.h
+)
 
 if (USE_PARALLEL)
     list (APPEND sm ${sm_parallel})

src/sm/Elements/LatticeElements/latticeframe3d.C

@@ -101,29 +101,29 @@ LatticeFrame3d::computeBmatrixAt(GaussPoint *aGaussPoint, FloatMatrix &answer, i
     answer.at(2, 3) =  0.;
     answer.at(2, 4) = 0.;
     answer.at(2, 5) = 0;
-    answer.at(2, 6) = -this->length*(1-this->s)/2.;
+    answer.at(2, 6) = -this->length*(1.-this->s)/2.;
     //Second node
     answer.at(2, 7) = 0.;
     answer.at(2, 8) = 1.;
     answer.at(2, 9) =  0.;
     answer.at(2, 10) = 0.;
     answer.at(2, 11) = 0;
-    answer.at(2, 12) = -this->length*(1+this->s)/2.;
+    answer.at(2, 12) = -this->length*(1.+this->s)/2.;
 
     //Shear displacement jump in z-plane
     //first node
     answer.at(3, 1) = 0.;
     answer.at(3, 2) = 0.;
     answer.at(3, 3) = -1.;
     answer.at(3, 4) = 0.;
-    answer.at(3, 5) = this->length*(1-this->s)/2.;
+    answer.at(3, 5) = this->length*(1.-this->s)/2.;
     answer.at(3, 6) = 0.;
     //Second node
     answer.at(3, 7) = 0.;
     answer.at(3, 8) = 0.;
     answer.at(3, 9) =  1.;
     answer.at(3, 10) = 0.;
-    answer.at(3, 11) = this->length*(1+this->s)/2.;
+    answer.at(3, 11) = this->length*(1.+this->s)/2.;
     answer.at(3, 12) = 0.;
 
     //Rotation around x-axis
@@ -241,7 +241,8 @@ LatticeFrame3d::computeStiffnessMatrix(FloatMatrix &answer, MatResponseMode rMod
     dbj.times(1. / length);
     bjt.beTranspositionOf(bj);
     answer.beProductOf(bjt, dbj);
-
+    //printf("answer/n");
+    //answer.printYourself();
     return;
 }
 
@@ -305,6 +306,7 @@ LatticeFrame3d::giveInternalForcesVector(FloatArray &answer,
     answer.clear();
 
     this->computeBmatrixAt(integrationRulesArray [ 0 ]->getIntegrationPoint(0), b);
+
     bt.beTranspositionOf(b);
 
     if ( useUpdatedGpRecord == 1 ) {
@@ -318,17 +320,13 @@ LatticeFrame3d::giveInternalForcesVector(FloatArray &answer,
 	strain.times(1./this->length);
         this->computeStressVector(stress, strain, integrationRulesArray [ 0 ]->getIntegrationPoint(0), tStep);
     }
-
     answer.beProductOf(bt, stress);
-
     if ( !this->isActivated(tStep) ) {
         answer.zero();
         return;
     }
 }
 
-
-
 bool
 LatticeFrame3d::computeGtoLRotationMatrix(FloatMatrix &answer)
 {
@@ -357,10 +355,10 @@ LatticeFrame3d::giveLocalCoordinateSystem(FloatMatrix &answer)
     Node *nodeA, *nodeB;
     nodeA = this->giveNode(1);
     nodeB = this->giveNode(2);
-
     lx.beDifferenceOf(nodeB->giveCoordinates(), nodeA->giveCoordinates() );
     lx.normalize();
 
+
     if ( this->referenceNode ) {
         Node *refNode = this->giveDomain()->giveNode(this->referenceNode);
         help.beDifferenceOf(refNode->giveCoordinates(), nodeA->giveCoordinates() );
@@ -446,6 +444,7 @@ LatticeFrame3d::initializeFrom(InputRecord &ir)
 
     this->s = 0.;
     IR_GIVE_OPTIONAL_FIELD(ir, s, _IFT_LatticeFrame3d_s);
+
 }
 
 
@@ -467,8 +466,6 @@ LatticeFrame3d::computeLength()
     return length;
 }
 
-
-
 void
 LatticeFrame3d::computeLumpedMassMatrix(FloatMatrix &answer, TimeStep *tStep)
 // Returns the lumped mass matrix of the receiver. This expression is

src/sm/Elements/LatticeElements/latticeframe3d.h

@@ -110,7 +110,7 @@ class LatticeFrame3d : public LatticeStructuralElement
 
 
 protected:
-    void computeBmatrixAt(GaussPoint *, FloatMatrix &, int = 1, int = ALL_STRAINS) override;
+    virtual void computeBmatrixAt(GaussPoint *, FloatMatrix &, int = 1, int = ALL_STRAINS) override;
     bool computeGtoLRotationMatrix(FloatMatrix &) override;
     void computeLumpedMassMatrix(FloatMatrix &answer, TimeStep *tStep) override;
     void computeMassMatrix(FloatMatrix &answer, TimeStep *tStep) override