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.

Sign up for GitHub

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

Fixed How SoilDyn is Called in the Glue Code #2

Merged
merged 1 commit into from
Apr 20, 2020
Merged

Fixed How SoilDyn is Called in the Glue Code #2

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
55 changes: 39 additions & 16 deletions docs/OtherSupporting/OpenFAST_Algorithms/OpenFAST_Algorithms.tex
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -44,6 +44,7 @@ \section{Definitions and Nomenclature}
AeroDyn & AD & AD \\
ServoDyn & SrvD & SrvD \\
SubDyn & SD & SD \\
ExtPtfm & ExtPtfm & ExtPtfm \\
HydroDyn & HydroDyn & HD \\
MAP++ & MAPp & MAP \\
FEAMooring & FEAM & FEAM \\
Expand Down Expand Up @@ -86,10 +87,10 @@ \section{Input-Output Relationships}
% SolveOption2a_Inp2BD
\State $\mathit{y\_ED} \gets \Call{ED\_CalcOutput}{\mathit{p\_ED},\mathit{u\_ED},\mathit{x\_ED},\mathit{xd\_ED},\mathit{z\_ED}}$
\State $\mathit{u\_BD} \gets \Call{TransferOutputsToInputs}{\mathit{y\_ED}}$
\State $\mathit{y\_BD} \gets \Call{BD\_CalcOutput}{\mathit{p\_BD},\mathit{u\_BD},\mathit{x\_BD},\mathit{xd\_BD},\mathit{z\_BD}}$

\State
% SolveOption2b_Inp2IfW
\State $\mathit{y\_BD} \gets \Call{BD\_CalcOutput}{\mathit{p\_BD},\mathit{u\_BD},\mathit{x\_BD},\mathit{xd\_BD},\mathit{z\_BD}}$
\State $\mathit{u\_AD}($no IfW$) \gets \Call{TransferOutputsToInputs}{\mathit{y\_ED}}$
\State $\mathit{u\_IfW} \gets \Call{TransferOutputsToInputs}{\mathit{y\_ED} at \mathit{u\_AD} nodes}$

Expand Down Expand Up @@ -120,11 +121,14 @@ \section{Input-Output Relationships}
% \State $\mathit{u\_BD} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$ % only if not BD_Solve_Option1
\State $\mathit{u\_HD} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_SD} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_ExtPtfm} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_MAP} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_FEAM} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_MD} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_Orca} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$

\State
\State $\mathit{y\_SD} \gets \Call{SD\_CalcOutput}{\mathit{p\_SD},\mathit{u\_SD},\mathit{x\_SD},\mathit{xd\_SD},\mathit{z\_SD}}$
\State $\mathit{u\_SlD} \gets \Call{TransferMeshPosition}{\mathit{y\_SD}}$
% \end Transfer_ED_to_HD_SD_BD_Mooring
%%%%
Expand All @@ -151,11 +155,11 @@ \section{Input-Output Relationships}


%\pagebreak %break here for now so that it doesn't look so strange
\subsection {Input-Output Solve for \textit{HydroDyn}, \textit{SubDyn}, \textit{MAP},
\textit{FEAMooring}, \textit{IceFloe}, and the Platform Reference Point Mesh in \textit{ElastoDyn}}
\subsection {Input-Output Solve for \textit{HydroDyn}, \textit{SubDyn}, \textit{OrcaFlexInterface}, \textit{BeamDyn}, \textit{SoilDyn}, \textit{ExtPtfm}, \textit{MAP}, \textit{FEAMooring}, \textit{MoorDyn},
\textit{FEAMooring}, \textit{IceFloe}, \textit{IceDyn}, and the Platform Reference Point Mesh in \textit{ElastoDyn}}

This procedure implements Solve Option 1 for the accelerations and loads in
\emph{HydroDyn},\emph{SubDyn},\emph{MAP},\emph{FEAMooring},\emph{OrcaFlexInterface},\emph{MoorDyn},\emph{SoilDyn}, and\emph{ElastoDyn} (at its platform reference point mesh).
\emph{HydroDyn},\emph{SubDyn},\emph{MAP},\emph{FEAMooring},\emph{OrcaFlexInterface},\emph{MoorDyn},\emph{SoilDyn}, \emph{BeamDyn}, \emph{ExtPtfm}, \emph{IceFloe}, \emph{IceDyn}, and \emph{ElastoDyn} (at its platform reference point mesh).
The other input-output relationships for these modules are solved using Solve Option 2.

%\begin{algorithm}[ht]
Expand All @@ -172,17 +176,17 @@ \section{Input-Output Relationships}
\State $\mathit{y\_IceD(:)} \gets \Call{CalcOutput}{\mathit{p\_IceD(:)},\mathit{u\_IceD(:)},\mathit{x\_IceD(:)},\mathit{xd\_IceD(:)},\mathit{z\_IceD(:)}}$
\State $\mathit{y\_SlD} \gets \Call{CalcOutput}{\mathit{p\_SlD},\mathit{u\_SlD},\mathit{x\_SlD},\mathit{xd\_SlD},\mathit{z\_SlD}}$
\State
\State\Comment{Form $u$ vector using loads and accelerations from $\mathit{u\_HD}$, $\mathit{u\_BD}$ $\mathit{u\_SD}$, and platform reference input from $\mathit{u\_ED}$}
\State\Comment{Form $u$ vector using loads and accelerations from $\mathit{u\_HD}$, $\mathit{u\_BD}$, $\mathit{u\_SD}$, \mathit{u\_Orca}, \mathit{u\_ExtPtfm}, and platform reference input from $\mathit{u\_ED}$}
\State
\State $u \gets \Call{u\_vec}{\mathit{u\_HD},\mathit{u\_SD},\mathit{u\_ED}}$
\State $u \gets \Call{u\_vec}{\mathit{u\_HD},\mathit{u\_SD},\mathit{u\_ED},\mathit{u\_BD},\mathit{u\_Orca},\mathit{u\_ExtPtfm}}$
\State $k \gets 0$
\Loop\Comment{Solve for loads and accelerations (direct feed-through terms)}
\State $y\_ED \gets \Call{ED\_CalcOutput}{\mathit{p\_ED},\mathit{u\_ED},\mathit{x\_ED},\mathit{xd\_ED},\mathit{z\_ED}}$
\State $y\_SD \gets \Call{SD\_CalcOutput}{\mathit{p\_SD},\mathit{u\_SD},\mathit{x\_SD},\mathit{xd\_SD},\mathit{z\_SD}}$
\State $y\_HD \gets \Call{HD\_CalcOutput}{\mathit{p\_HD},\mathit{u\_HD},\mathit{x\_HD},\mathit{xd\_HD},\mathit{z\_HD}}$
\State $y\_BD \gets \Call{BD\_CalcOutput}{\mathit{p\_BD},\mathit{u\_BD},\mathit{x\_BD},\mathit{xd\_BD},\mathit{z\_BD}}$
\State $y\_Orca \gets \Call{Orca\_CalcOutput}{\mathit{p\_Orca},\mathit{u\_Orca},\mathit{x\_Orca},\mathit{xd\_Orca},\mathit{z\_Orca}}$
\State $y\_SlD \gets \Call{SlD\_CalcOutput}{\mathit{p\_SlD},\mathit{u\_SlD},\mathit{x\_SlD},\mathit{xd\_SlD},\mathit{z\_SlD}}$
\State $\mathit{y\_ExtPtfm} \gets \Call{CalcOutput}{\mathit{p\_ExtPtfm},\mathit{u\_ExtPtfm},\mathit{x\_ExtPtfm},\mathit{xd\_ExtPtfm},\mathit{z\_ExtPtfm}}$

\If{ $k \geq k\_max$}
\State exit loop
Expand All @@ -191,6 +195,8 @@ \section{Input-Output Relationships}
\State$\mathit{u\_BD\_tmp} \gets \Call{TransferMeshMotions}{y\_ED}$
\State$\mathit{u\_MAP\_tmp} \gets \Call{TransferMeshMotions}{y\_ED}$
\State$\mathit{u\_FEAM\_tmp} \gets \Call{TransferMeshMotions}{y\_ED}$
\State$\mathit{u\_Orca\_tmp} \gets \Call{TransferMeshMotions}{y\_ED}$
\State$\mathit{u\_MD\_tmp} \gets \Call{TransferMeshMotions}{y\_ED}$
\State$\mathit{u\_IceF\_tmp} \gets \Call{TransferMeshMotions}{y\_SD}$
\State$\mathit{u\_IceD\_tmp(:)} \gets \Call{TransferMeshMotions}{y\_SD}$
\State$\mathit{u\_SlD\_tmp} \gets \Call{TransferMeshMotions}{y\_SD}$
Expand All @@ -216,7 +222,7 @@ \section{Input-Output Relationships}
\end{aligned}$

\State
\State$\mathit{U\_Residual} \gets u - \Call{u\_vec}{\mathit{u\_HD\_tmp},\mathit{u\_SD\_tmp},\mathit{u\_ED\_tmp},\mathit{u\_SlD\_tmp}}$
\State$\mathit{U\_Residual} \gets u - \Call{u\_vec}{\mathit{u\_HD\_tmp},\mathit{u\_SD\_tmp},\mathit{u\_ED\_tmp},\mathit{u\_BD\_tmp},\mathit{u\_Orca\_tmp},\mathit{u\_ExtPtfm\_tmp}}$
\State

\If{ last Jacobian was calculated at least $\mathit{DT\_UJac}$ seconds ago }
Expand All @@ -239,16 +245,18 @@ \section{Input-Output Relationships}
\EndIf
\State
\State $u \gets u + \Delta u$
\State Transfer $u$ to $\mathit{u\_HD}$, $\mathit{u\_SD}$, and $\mathit{u\_ED}$\Comment{loads and accelerations only}
\State Transfer $u$ to $\mathit{u\_HD}$, $\mathit{u\_SD}$, $\mathit{u\_BD}$, $\mathit{u\_Orca}$, $\mathit{u\_ExtPtfm}$, and $\mathit{u\_ED}$\Comment{loads and accelerations only}
\State $k=k+1$

\EndLoop

\State\Comment{Transfer non-acceleration fields to motion input meshes}
\State

\State$\mathit{u\_BD}($not accelerations$) \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State$\mathit{u\_HD}($not accelerations$) \gets \Call{TransferMeshMotions}{\mathit{y\_ED},\mathit{y\_SD}}$
\State$\mathit{u\_SD}($not accelerations$) \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State$\mathit{u\_Orca}($not accelerations$) \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$ \State$\mathit{u\_ExtPtfm}($not accelerations$) \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State
\State $\mathit{u\_MAP} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_MD} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
Expand All @@ -260,13 +268,6 @@ \section{Input-Output Relationships}
\EndProcedure
\end{algorithmic}

\paragraph{SoilDyn}
The motion inputs to \emph{SoilDyn} are the position and orientation of the mesh points. The \emph{SoilDyn} module expects changes in position to be in the order of $10^{-3}$~m, or angular changes of $10^{-4}$~radians or less throughout the duration of the simulation with response terms on the order of $10^{12}$~N/m (or N-m/radians).
To simplify the mathematics in the Jacobian, an Euler angle extraction is performed to get rotations about the $x$, $y$, and $z$ axis (the angles are small enough that order does not matter in this context). So the $u$ input vector for \emph{SoilDyn} is a 6-vector for each mesh point.
The resulting reaction forces from \emph{SoilDyn} are roughly characterized as $F_{\mathrm{react}} \approx k x$ where $k$ is a stiffness matrix and $x$ is the 6-vector of translational and rotational displacements (note that there may be some damping or drift terms, but of much lower magnitude).
This creates a situation where $\mathit{u_SD} \gets \mathit{y_SlD} \gets \mathit{y_SD}$ is a feed back loop of \emph{SubDyn} node displacements (stored in states) of a given mesh point directly through \emph{SoilDyn} to the input force on that node point. Without including this calculation in the Jacobian, the accelerations of the \emph{SubDyn} nodes will be missing the reaction force and likely cause non-convergence or large oscillations. For this reason we include the \emph{SoilDyn} input / output into the Jacobian. However this does cause some of the Jacobian to differ by as many as six orders of magnitude, which is numerically non-ideal for the solve.



\subsection {Implementation of line2-to-line2 loads mapping}
The inverse-lumping of loads is computed by a block matrix solve for the distributed forces and moments,
Expand Down Expand Up @@ -351,6 +352,28 @@ \section{Solve Option 2 Improvements}
\State $\Call{AD\_UpdateStates}{\mathit{p\_AD},\mathit{u\_AD},\mathit{x\_AD},\mathit{xd\_AD},\mathit{z\_AD}}$
\State $\Call{SrvD\_UpdateStates}{\mathit{p\_SrvD},\mathit{u\_SrvD},\mathit{x\_SrvD},\mathit{xd\_SrvD},\mathit{z\_SrvD}}$
\State
% \begin Transfer_ED_to_HD_SD_BD_Mooring
% \State $\mathit{u\_ED}($not platform reference point$) \gets \Call{TransferOutputsToInputs}{y\_SrvD,y\_AD}$ %\Comment{sets all but platform reference point inputs}

% \State $\mathit{u\_BD} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$ % only if not BD_Solve_Option1
\State $\mathit{u\_HD} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_SD} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_ExtPtfm} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_MAP} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_FEAM} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_MD} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
\State $\mathit{u\_Orca} \gets \Call{TransferMeshMotions}{\mathit{y\_ED}}$
% \end Transfer_ED_to_HD_SD_BD_Mooring
\State
\State $\Call{SD\_UpdateStates}{\mathit{p\_SD},\mathit{u\_SD},\mathit{x\_SD},\mathit{xd\_SD},\mathit{z\_SD}}$
\State $\mathit{y\_SD} \gets \Call{SD\_CalcOutput}{\mathit{p\_SD},\mathit{u\_SD},\mathit{x\_SD},\mathit{xd\_SD},\mathit{z\_SD}}$
\State $\mathit{u\_SlD} \gets \Call{TransferMeshPosition}{\mathit{y\_SD}}$
\State $\Call{SlD\_UpdateStates}{\mathit{p\_SlD},\mathit{u\_SlD},\mathit{x\_SlD},\mathit{xd\_SlD},\mathit{z\_SlD}}$

%%%%

\State

\State All other modules (used in Solve Option 1) advance their states
\EndProcedure
\end{algorithmic}
Expand Down
2 changes: 1 addition & 1 deletion modules/openfast-library/src/FAST_Solver.f90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5088,7 +5088,7 @@ SUBROUTINE SolveOption2b_Inp2IfW(this_time, this_state, p_FAST, m_FAST, ED, BD,
ErrMsg = ""


IF ( p_FAST%CompElast == Module_BD .AND. .NOT. BD_Solve_Option1 ) THEN
IF ( p_FAST%CompElast == Module_BD ) THEN
DO k=1,p_FAST%nBeams
CALL BD_CalcOutput( this_time, BD%Input(1,k), BD%p(k), BD%x(k,this_state), BD%xd(k,this_state),&
BD%z(k,this_state), BD%OtherSt(k,this_state), BD%y(k), BD%m(k), ErrStat2, ErrMsg2 )
Expand Down