![]() |
OpenFAST
Wind turbine multiphysics simulator
|
This module contains the routines used by FAST to solve input-output equations and to advance states. More...
Functions/Subroutines | |
subroutine | bd_inputsolve (p_FAST, BD, y_AD, u_AD, y_ED, MeshMapData, ErrStat, ErrMsg) |
This routine sets the inputs required for BD–using the Option 2 solve method; currently the only inputs solved in this routine are the blade distributed loads from AD15; other inputs are solved in option 1. More... | |
subroutine | ed_inputsolve (p_FAST, u_ED, y_ED, p_AD14, y_AD14, y_AD, y_SrvD, u_AD, u_SrvD, MeshMapData, ErrStat, ErrMsg) |
This routine sets the inputs required for ED–using the Option 2 solve method; currently the only input not solved in this routine are the fields on PlatformPtMesh and HubPtLoad, which are solved in option 1. More... | |
subroutine | ifw_inputsolve (p_FAST, m_FAST, u_IfW, p_IfW, u_AD14, u_AD, OtherSt_AD, y_ED, ErrStat, ErrMsg) |
This routine determines the points in space where InflowWind needs to compute wind speeds. More... | |
subroutine | ifw_setexternalinputs (p_IfW, m_FAST, y_ED, u_IfW) |
This routine sets the inputs required for InflowWind from an external source (Simulink) More... | |
subroutine | ad_inputsolve_ifw (p_FAST, u_AD, y_IfW, y_OpFM, ErrStat, ErrMsg) |
This routine sets the AeroDyn wind inflow inputs. More... | |
subroutine | ad_inputsolve_noifw (p_FAST, u_AD, y_SrvD, y_ED, BD, MeshMapData, ErrStat, ErrMsg) |
This routine sets all the AeroDyn inputs, except for the wind inflow values. More... | |
subroutine | ad14_inputsolve_ifw (p_FAST, u_AD14, y_IfW, y_OpFM, ErrStat, ErrMsg) |
This routine sets the AeroDyn14 wind inflow inputs. More... | |
subroutine | ad14_inputsolve_noifw (p_FAST, u_AD14, y_ED, MeshMapData, ErrStat, ErrMsg) |
This routine sets all the AeroDyn14 inputs, except for the wind inflow values. More... | |
subroutine | srvd_inputsolve (p_FAST, m_FAST, u_SrvD, y_ED, y_IfW, y_OpFM, y_BD, MeshMapData, ErrStat, ErrMsg) |
This routine sets the inputs required for ServoDyn. More... | |
subroutine | srvd_setexternalinputs (p_FAST, m_FAST, u_SrvD) |
This routine sets the inputs required for ServoDyn from an external source (Simulink) More... | |
subroutine | transfer_sd_to_hd (y_SD, u_HD_M_LumpedMesh, u_HD_M_DistribMesh, MeshMapData, ErrStat, ErrMsg) |
This routine transfers the SD outputs into inputs required for HD. More... | |
subroutine | transfer_platformmotion_to_hd (PlatformMotion, u_HD, MeshMapData, ErrStat, ErrMsg) |
This routine transfers the platform motion output of the structural module (ED) into inputs required for HD. More... | |
subroutine | transfer_ed_to_hd_sd_bd_mooring (p_FAST, y_ED, u_HD, u_SD, u_ExtPtfm, u_MAP, u_FEAM, u_MD, u_Orca, u_BD, MeshMapData, ErrStat, ErrMsg) |
This routine transfers the ED outputs into inputs required for HD, SD, ExtPtfm, BD, MAP, and/or FEAM. More... | |
subroutine | map_inputsolve (u_MAP, y_ED, MeshMapData, ErrStat, ErrMsg) |
This routine sets the inputs required for MAP. More... | |
subroutine | feam_inputsolve (u_FEAM, y_ED, MeshMapData, ErrStat, ErrMsg) |
This routine sets the inputs required for FEAM. More... | |
subroutine | md_inputsolve (u_MD, y_ED, MeshMapData, ErrStat, ErrMsg) |
This routine sets the inputs required for MoorDyn. More... | |
subroutine | icefloe_inputsolve (u_IceF, y_SD, MeshMapData, ErrStat, ErrMsg) |
This routine sets the inputs required for IceFloe. More... | |
subroutine | iced_inputsolve (u_IceD, y_SD, MeshMapData, legNum, ErrStat, ErrMsg) |
This routine sets the inputs required for IceFloe. More... | |
subroutine | transfer_ed_to_bd (y_ED, u_BD, MeshMapData, ErrStat, ErrMsg) |
This routine sets the inputs required for IceFloe. More... | |
subroutine | transfer_ed_to_bd_tmp (y_ED, MeshMapData, ErrStat, ErrMsg) |
This routine sets the inputs required for IceFloe. More... | |
subroutine | transfer_hd_to_sd (u_mapped, u_SD_LMesh, u_mapped_positions, y_HD, u_HD_M_LumpedMesh, u_HD_M_DistribMesh, MeshMapData, ErrStat, ErrMsg) |
This routine transfers the HD outputs into inputs required for ED. More... | |
real(reki) function | getperturb (x) |
function to return the size of perturbation in calculating jacobian with finite differences. More... | |
subroutine | ed_hd_inputoutputsolve (this_time, p_FAST, calcJacobian, u_ED, p_ED, x_ED, xd_ED, z_ED, OtherSt_ED, y_ED, m_ED, u_HD, p_HD, x_HD, xd_HD, z_HD, OtherSt_HD, y_HD, m_HD, u_MAP, y_MAP, u_FEAM, y_FEAM, u_MD, y_MD, MeshMapData, ErrStat, ErrMsg, WriteThisStep) |
This routine performs the Input-Output solve for ED and HD. More... | |
subroutine | fullopt1_inputoutputsolve (this_time, p_FAST, calcJacobian, u_ED, p_ED, x_ED, xd_ED, z_ED, OtherSt_ED, y_ED, m_ED, u_SD, p_SD, x_SD, xd_SD, z_SD, OtherSt_SD, y_SD, m_SD, u_ExtPtfm, p_ExtPtfm, x_ExtPtfm, xd_ExtPtfm, z_ExtPtfm, OtherSt_ExtPtfm, y_ExtPtfm, m_ExtPtfm, u_HD, p_HD, x_HD, xd_HD, z_HD, OtherSt_HD, y_HD, m_HD, u_BD, p_BD, x_BD, xd_BD, z_BD, OtherSt_BD, y_BD, m_BD, u_Orca, p_Orca, x_Orca, xd_Orca, z_Orca, OtherSt_Orca, y_Orca, m_Orca, u_MAP, y_MAP, u_FEAM, y_FEAM, u_MD, y_MD, u_IceF, y_IceF, u_IceD, y_IceD, MeshMapData, ErrStat, ErrMsg, WriteThisStep) |
This routine performs the Input-Output solve for ED, SD, HD, BD, and/or the OrcaFlex Interface. More... | |
subroutine | init_fullopt1_jacobian (p_FAST, MeshMapData, ED_PlatformPtMesh, SD_TPMesh, SD_LMesh, HD_M_LumpedMesh, HD_M_DistribMesh, HD_WAMIT_Mesh, ED_HubPtLoad, u_BD, Orca_PtfmMesh, ExtPtfm_PtfmMesh, ErrStat, ErrMsg) |
This routine initializes the array that maps rows/columns of the Jacobian to specific mesh fields. More... | |
subroutine | create_fullopt1_uvector (u, ED_PlatformPtMesh, SD_TPMesh, SD_LMesh, HD_M_LumpedMesh, HD_M_DistribMesh, HD_WAMIT_Mesh, ED_HubPtLoad, BD_RootMotion, Orca_PtfmMesh, ExtPtfm_PtfmMesh, p_FAST) |
This routine basically packs the relevant parts of the modules' input meshes for use in this InputOutput solve. More... | |
subroutine | add_fullopt1_u_delta (p_FAST, Jac_u_indx, u_delta, u_ED, u_SD, u_HD, u_BD, u_Orca, u_ExtPtfm) |
This routine adds u_delta to the corresponding mesh field and scales it as appropriate. More... | |
subroutine | perturb_u_fullopt1 (p_FAST, Jac_u_indx, n, u_perturb, u_ED_perturb, u_SD_perturb, u_HD_perturb, u_BD_perturb, u_Orca_perturb, u_ExtPtfm_perturb, perturb) |
This routine perturbs the nth element of the u array (and mesh/field it corresponds to) More... | |
subroutine | resetremapflags (p_FAST, ED, BD, AD14, AD, HD, SD, ExtPtfm, SrvD, MAPp, FEAM, MD, Orca, IceF, IceD) |
This routine resets the remap flags on all of the meshes. More... | |
subroutine | initmodulemappings (p_FAST, ED, BD, AD14, AD, HD, SD, ExtPtfm, SrvD, MAPp, FEAM, MD, Orca, IceF, IceD, MeshMapData, ErrStat, ErrMsg) |
This routine initializes all of the mapping data structures needed between the various modules. More... | |
subroutine | calcoutputs_and_solveforinputs (n_t_global, this_time, this_state, calcJacobian, NextJacCalcTime, p_FAST, m_FAST, WriteThisStep, ED, BD, SrvD, AD14, AD, IfW, OpFM, HD, SD, ExtPtfm, MAPp, FEAM, MD, Orca, IceF, IceD, MeshMapData, ErrStat, ErrMsg) |
This subroutine solves the input-output relations for all of the modules. More... | |
subroutine | solveoption1 (this_time, this_state, calcJacobian, p_FAST, ED, BD, HD, SD, ExtPtfm, MAPp, FEAM, MD, Orca, IceF, IceD, MeshMapData, ErrStat, ErrMsg, WriteThisStep) |
This routine implements the "option 1" solve for all inputs with direct links to HD, SD, ExtPtfm, MAP, OrcaFlex interface, and the ED platform reference point. More... | |
subroutine | solveoption2a_inp2bd (this_time, this_state, p_FAST, m_FAST, ED, BD, AD14, AD, SrvD, IfW, OpFM, MeshMapData, ErrStat, ErrMsg, WriteThisStep) |
This routine implements the first part of the "option 2" solve for inputs that apply to BeamDyn and AeroDyn. More... | |
subroutine | solveoption2b_inp2ifw (this_time, this_state, p_FAST, m_FAST, ED, BD, AD14, AD, SrvD, IfW, OpFM, MeshMapData, ErrStat, ErrMsg, WriteThisStep) |
This routine implements the first part of the "option 2" solve for inputs that apply to AeroDyn & InflowWind. More... | |
subroutine | solveoption2c_inp2ad_srvd (this_time, this_state, p_FAST, m_FAST, ED, BD, AD14, AD, SrvD, IfW, OpFM, MeshMapData, ErrStat, ErrMsg, WriteThisStep) |
This routine implements the first part of the "option 2" solve for inputs that apply to AeroDyn and ServoDyn. More... | |
subroutine | solveoption2 (this_time, this_state, p_FAST, m_FAST, ED, BD, AD14, AD, SrvD, IfW, OpFM, MeshMapData, ErrStat, ErrMsg, firstCall, WriteThisStep) |
This routine implements the "option 2" solve for all inputs without direct links to HD, SD, MAP, or the ED platform reference point. More... | |
subroutine | fast_advancestates (t_initial, n_t_global, p_FAST, m_FAST, ED, BD, SrvD, AD14, AD, IfW, OpFM, HD, SD, ExtPtfm, MAPp, FEAM, MD, Orca, IceF, IceD, MeshMapData, ErrStat, ErrMsg, WriteThisStep) |
This routines advances the states of each module. More... | |
subroutine | fast_extrapinterpmods (t_global_next, p_FAST, m_FAST, ED, BD, SrvD, AD14, AD, IfW, HD, SD, ExtPtfm, MAPp, FEAM, MD, Orca, IceF, IceD, ErrStat, ErrMsg) |
This routine extrapolates inputs to modules to give predicted values at t+dt. More... | |
This module contains the routines used by FAST to solve input-output equations and to advance states.
subroutine fast_solver::ad14_inputsolve_ifw | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(ad14_inputtype), intent(inout) | u_AD14, | ||
type(inflowwind_outputtype), intent(in) | y_IfW, | ||
type(opfm_outputtype), intent(in) | y_OpFM, | ||
integer(intki) | ErrStat, | ||
character(*) | ErrMsg | ||
) |
This routine sets the AeroDyn14 wind inflow inputs.
[in] | p_fast | parameter FAST data |
[in,out] | u_ad14 | The inputs to AeroDyn14 |
[in] | y_ifw | The outputs from InflowWind |
[in] | y_opfm | outputs from the OpenFOAM integration module |
errstat | Error status of the operation | |
errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::ad14_inputsolve_noifw | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(ad14_inputtype), intent(inout) | u_AD14, | ||
type(ed_outputtype), intent(in) | y_ED, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki) | ErrStat, | ||
character(*) | ErrMsg | ||
) |
This routine sets all the AeroDyn14 inputs, except for the wind inflow values.
THIS ROUTINE IS A HACK TO GET THE OUTPUTS FROM ELASTODYN INTO AERODYN14. DO NOT COPY OR USE IN NEW CODE!
[in] | p_fast | parameter FAST data |
[in,out] | u_ad14 | The inputs to AeroDyn14 |
[in] | y_ed | The outputs from the structural dynamics module |
[in,out] | meshmapdata | Data for mapping between modules |
errstat | Error status of the operation | |
errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::ad_inputsolve_ifw | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(ad_inputtype), intent(inout) | u_AD, | ||
type(inflowwind_outputtype), intent(in) | y_IfW, | ||
type(opfm_outputtype), intent(in) | y_OpFM, | ||
integer(intki) | ErrStat, | ||
character(*) | ErrMsg | ||
) |
This routine sets the AeroDyn wind inflow inputs.
[in] | p_fast | FAST parameter data |
[in,out] | u_ad | The inputs to AeroDyn |
[in] | y_ifw | The outputs from InflowWind |
[in] | y_opfm | outputs from the OpenFOAM integration module |
errstat | Error status of the operation | |
errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::ad_inputsolve_noifw | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(ad_inputtype), intent(inout) | u_AD, | ||
type(srvd_outputtype), intent(in) | y_SrvD, | ||
type(ed_outputtype), intent(in) | y_ED, | ||
type(beamdyn_data), intent(in) | BD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki) | ErrStat, | ||
character(*) | ErrMsg | ||
) |
This routine sets all the AeroDyn inputs, except for the wind inflow values.
[in] | p_fast | FAST parameter data |
[in,out] | u_ad | The inputs to AeroDyn14 |
[in] | y_srvd | ServoDyn outputs |
[in] | y_ed | The outputs from the structural dynamics module |
[in] | bd | The data from BeamDyn (want the outputs only, but it's in an array) |
[in,out] | meshmapdata | Data for mapping between modules |
errstat | Error status of the operation | |
errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::add_fullopt1_u_delta | ( | type(fast_parametertype), intent(in) | p_FAST, |
integer( intki ), dimension(:,:), intent(in) | Jac_u_indx, | ||
real( reki ), dimension(:), intent(in) | u_delta, | ||
type(ed_inputtype), intent(inout) | u_ED, | ||
type(sd_inputtype), intent(inout) | u_SD, | ||
type(hydrodyn_inputtype), intent(inout) | u_HD, | ||
type(bd_inputtype), dimension(:), intent(inout) | u_BD, | ||
type(orca_inputtype), intent(inout) | u_Orca, | ||
type(extptfm_inputtype), intent(inout) | u_ExtPtfm | ||
) |
This routine adds u_delta to the corresponding mesh field and scales it as appropriate.
[in] | p_fast | Glue-code simulation parameters |
[in] | jac_u_indx | Index to map Jacobian u-vector into mesh fields |
[in] | u_delta | The delta amount to add to the appropriate mesh fields |
[in,out] | u_ed | ED System inputs |
[in,out] | u_sd | SD System inputs |
[in,out] | u_hd | SD System inputs |
[in,out] | u_bd | BD System inputs |
[in,out] | u_orca | Orca System inputs |
[in,out] | u_extptfm | ExtPtfm System inputs |
subroutine fast_solver::bd_inputsolve | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(beamdyn_data), intent(inout) | BD, | ||
type(ad_outputtype), intent(in) | y_AD, | ||
type(ad_inputtype), intent(in) | u_AD, | ||
type(ed_outputtype), intent(in) | y_ED, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine sets the inputs required for BD–using the Option 2 solve method; currently the only inputs solved in this routine are the blade distributed loads from AD15; other inputs are solved in option 1.
[in] | p_fast | Glue-code simulation parameters |
[in,out] | bd | BD Inputs at t |
[in] | y_ad | AeroDyn outputs |
[in] | u_ad | AD inputs (for AD-BD load transfer) |
[in] | y_ed | ElastoDyn outputs |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status |
[out] | errmsg | Error message |
subroutine fast_solver::calcoutputs_and_solveforinputs | ( | integer(intki), intent(in) | n_t_global, |
real(dbki), intent(in) | this_time, | ||
integer(intki), intent(in) | this_state, | ||
logical, intent(inout) | calcJacobian, | ||
real(dbki), intent(in) | NextJacCalcTime, | ||
type(fast_parametertype), intent(in) | p_FAST, | ||
type(fast_miscvartype), intent(in) | m_FAST, | ||
logical, intent(in) | WriteThisStep, | ||
type(elastodyn_data), intent(inout) | ED, | ||
type(beamdyn_data), intent(inout) | BD, | ||
type(servodyn_data), intent(inout) | SrvD, | ||
type(aerodyn14_data), intent(inout) | AD14, | ||
type(aerodyn_data), intent(inout) | AD, | ||
type(inflowwind_data), intent(inout) | IfW, | ||
type(openfoam_data), intent(inout) | OpFM, | ||
type(hydrodyn_data), intent(inout) | HD, | ||
type(subdyn_data), intent(inout) | SD, | ||
type(extptfm_data), intent(inout) | ExtPtfm, | ||
type(map_data), intent(inout) | MAPp, | ||
type(feamooring_data), intent(inout) | FEAM, | ||
type(moordyn_data), intent(inout) | MD, | ||
type(orcaflex_data), intent(inout) | Orca, | ||
type(icefloe_data), intent(inout) | IceF, | ||
type(icedyn_data), intent(inout) | IceD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This subroutine solves the input-output relations for all of the modules.
It is a subroutine because it gets done twice– once at the start of the n_t_global loop and once in the j_pc loop, using different states. *** Note that modules that do not have direct feedthrough should be called first. ***
[in] | this_time | The current simulation time (actual or time of prediction) |
[in] | this_state | Index into the state array (current or predicted states) |
[in] | n_t_global | current time step (used only for SrvD hack) |
[in,out] | calcjacobian | Should we calculate Jacobians in Option 1? |
[in] | nextjaccalctime | Time between calculating Jacobians in the HD-ED and SD-ED simulations |
[in] | p_fast | Parameters for the glue code |
[in] | m_fast | Misc variables (including external inputs) for the glue code |
[in] | writethisstep | Will we print the WriteOutput values this step? |
[in,out] | ed | ElastoDyn data |
[in,out] | bd | BeamDyn data |
[in,out] | srvd | ServoDyn data |
[in,out] | ad14 | AeroDyn14 data |
[in,out] | ad | AeroDyn data |
[in,out] | ifw | InflowWind data |
[in,out] | opfm | OpenFOAM data |
[in,out] | hd | HydroDyn data |
[in,out] | sd | SubDyn data |
[in,out] | extptfm | ExtPtfm data |
[in,out] | mapp | MAP data |
[in,out] | feam | FEAMooring data |
[in,out] | md | Data for the MoorDyn module |
[in,out] | orca | OrcaFlex interface data |
[in,out] | icef | IceFloe data |
[in,out] | iced | All the IceDyn data used in time-step loop |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
For cases with HydroDyn, BeamDyn, OrcaFlex interface, and/or SubDyn, it calls ED_CalcOuts (a time-sink) 3 times per step/correction (plus the 6 calls when calculating the Jacobian). In cases without HydroDyn or SubDyn, it is the same as Option 1 before 2 (with 1 call to ED_CalcOuts either way).
Option 1 before 2 usually requires a correction step, whereas Option 2 before Option 1 often does not. Thus we are using this option, calling ED_CalcOuts 3 times (option 2 before 1 with no correction step) instead of 4 times (option1 before 2 with one correction step). Note that this analysis may change if/when AeroDyn14 (and ServoDyn?) generate different outputs on correction steps. (Currently, AeroDyn (v14) returns old values until time advances.) Also note that AD15's DBEMT module (and UA?) is heavily time-dependent without calling the structural code first (DBEMT's filters do not deal well with the extrapolated inputs).
Solve option 2 (modules without direct feedthrough):
transfer ED outputs to other modules used in option 1:
Solve option 1 (rigorous solve on loads/accelerations)
Now use the ElastoDyn and BD outputs from option1 to update the inputs for InflowWind, AeroDyn, and ServoDyn (necessary only if they have states)
subroutine fast_solver::create_fullopt1_uvector | ( | real(reki), dimension(:), intent(inout) | u, |
type(meshtype), intent(in) | ED_PlatformPtMesh, | ||
type(meshtype), intent(in) | SD_TPMesh, | ||
type(meshtype), intent(in) | SD_LMesh, | ||
type(meshtype), intent(in) | HD_M_LumpedMesh, | ||
type(meshtype), intent(in) | HD_M_DistribMesh, | ||
type(meshtype), intent(in) | HD_WAMIT_Mesh, | ||
type(meshtype), intent(in) | ED_HubPtLoad, | ||
type(meshtype), dimension(:), intent(in) | BD_RootMotion, | ||
type(meshtype), intent(in) | Orca_PtfmMesh, | ||
type(meshtype), intent(in) | ExtPtfm_PtfmMesh, | ||
type(fast_parametertype), intent(in) | p_FAST | ||
) |
This routine basically packs the relevant parts of the modules' input meshes for use in this InputOutput solve.
Do not change the order of this packing without changing subroutine Init_FullOpt1_Jacobian()!
[in,out] | u | output u vector |
[in] | ed_platformptmesh | ElastoDyn PlatformPt mesh |
[in] | sd_tpmesh | SubDyn TP mesh |
[in] | sd_lmesh | SubDyn Lmesh |
[in] | hd_m_lumpedmesh | HydroDyn Morison Lumped mesh |
[in] | hd_m_distribmesh | HydroDyn Morison distributed mesh |
[in] | hd_wamit_mesh | HydroDyn WAMIT mesh |
[in] | ed_hubptload | ElastoDyn HubPt mesh |
[in] | bd_rootmotion | BeamDyn RootMotion meshes |
[in] | orca_ptfmmesh | OrcaFlex interface PtfmMesh |
[in] | extptfm_ptfmmesh | ExtPtfm interface PtfmMesh |
[in] | p_fast | FAST parameters |
subroutine fast_solver::ed_hd_inputoutputsolve | ( | real(dbki), intent(in) | this_time, |
type(fast_parametertype), intent(in) | p_FAST, | ||
logical, intent(in) | calcJacobian, | ||
type(ed_inputtype), intent(inout) | u_ED, | ||
type(ed_parametertype), intent(in) | p_ED, | ||
type(ed_continuousstatetype), intent(in) | x_ED, | ||
type(ed_discretestatetype), intent(in) | xd_ED, | ||
type(ed_constraintstatetype), intent(in) | z_ED, | ||
type(ed_otherstatetype), intent(inout) | OtherSt_ED, | ||
type(ed_outputtype), intent(inout) | y_ED, | ||
type(ed_miscvartype), intent(inout) | m_ED, | ||
type(hydrodyn_inputtype), intent(inout) | u_HD, | ||
type(hydrodyn_parametertype), intent(in) | p_HD, | ||
type(hydrodyn_continuousstatetype), intent(in) | x_HD, | ||
type(hydrodyn_discretestatetype), intent(in) | xd_HD, | ||
type(hydrodyn_constraintstatetype), intent(in) | z_HD, | ||
type(hydrodyn_otherstatetype), intent(inout) | OtherSt_HD, | ||
type(hydrodyn_outputtype), intent(inout) | y_HD, | ||
type(hydrodyn_miscvartype), intent(inout) | m_HD, | ||
type(map_inputtype), intent(inout) | u_MAP, | ||
type(map_outputtype), intent(in) | y_MAP, | ||
type(feam_inputtype), intent(inout) | u_FEAM, | ||
type(feam_outputtype), intent(in) | y_FEAM, | ||
type(md_inputtype), intent(inout) | u_MD, | ||
type(md_outputtype), intent(in) | y_MD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
logical, intent(in) | WriteThisStep | ||
) |
This routine performs the Input-Output solve for ED and HD.
Note that this has been customized for the physics in the problems and is not a general solution.
[in] | this_time | The current simulation time (actual or time of prediction) |
[in] | p_fast | Glue-code simulation parameters |
[in] | calcjacobian | Should we calculate Jacobians this time? (should be TRUE on initialization, then can be false [significantly reducing computational time]) |
[in] | x_ed | Continuous states |
[in] | xd_ed | Discrete states |
[in] | z_ed | Constraint states |
[in,out] | otherst_ed | Other states |
[in] | p_ed | Parameters |
[in,out] | u_ed | System inputs |
[in,out] | y_ed | System outputs |
[in,out] | m_ed | misc/optimization variables |
[in] | x_hd | Continuous states |
[in] | xd_hd | Discrete states |
[in] | z_hd | Constraint states |
[in,out] | otherst_hd | Other states |
[in] | p_hd | Parameters |
[in,out] | u_hd | System inputs |
[in,out] | y_hd | System outputs |
[in,out] | m_hd | misc/optimization variables |
[in] | y_map | MAP outputs |
[in,out] | u_map | MAP inputs (INOUT just because I don't want to use another tempoarary mesh and we'll overwrite this later) |
[in] | y_feam | FEAM outputs |
[in,out] | u_feam | FEAM inputs (INOUT just because I don't want to use another tempoarary mesh and we'll overwrite this later) |
[in] | y_md | MoorDyn outputs |
[in,out] | u_md | MoorDyn inputs (INOUT just because I don't want to use another tempoarary mesh and we'll overwrite this later) |
[in,out] | meshmapdata | data for mapping meshes between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
[in] | writethisstep | Will we print the WriteOutput values this step? |
subroutine fast_solver::ed_inputsolve | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(ed_inputtype), intent(inout) | u_ED, | ||
type(ed_outputtype), intent(in) | y_ED, | ||
type(ad14_parametertype), intent(in) | p_AD14, | ||
type(ad14_outputtype), intent(in) | y_AD14, | ||
type(ad_outputtype), intent(in) | y_AD, | ||
type(srvd_outputtype), intent(in) | y_SrvD, | ||
type(ad_inputtype), intent(in) | u_AD, | ||
type(srvd_inputtype), intent(in) | u_SrvD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine sets the inputs required for ED–using the Option 2 solve method; currently the only input not solved in this routine are the fields on PlatformPtMesh and HubPtLoad, which are solved in option 1.
[in] | p_fast | Glue-code simulation parameters |
[in,out] | u_ed | ED Inputs at t |
[in] | y_ed | ElastoDyn outputs (need translation displacement on meshes for loads mapping) |
[in] | p_ad14 | AeroDyn14 parameters (a hack because the AD14 meshes aren't set up properly) |
[in] | y_ad14 | AeroDyn14 outputs |
[in] | y_ad | AeroDyn outputs |
[in] | u_ad | AD inputs (for AD-ED load transfer) |
[in] | y_srvd | ServoDyn outputs |
[in] | u_srvd | ServoDyn inputs |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status |
[out] | errmsg | Error message |
subroutine fast_solver::fast_advancestates | ( | real(dbki), intent(in) | t_initial, |
integer(intki), intent(in) | n_t_global, | ||
type(fast_parametertype), intent(in) | p_FAST, | ||
type(fast_miscvartype), intent(in) | m_FAST, | ||
type(elastodyn_data), intent(inout) | ED, | ||
type(beamdyn_data), intent(inout) | BD, | ||
type(servodyn_data), intent(inout) | SrvD, | ||
type(aerodyn14_data), intent(inout) | AD14, | ||
type(aerodyn_data), intent(inout) | AD, | ||
type(inflowwind_data), intent(inout) | IfW, | ||
type(openfoam_data), intent(inout) | OpFM, | ||
type(hydrodyn_data), intent(inout) | HD, | ||
type(subdyn_data), intent(inout) | SD, | ||
type(extptfm_data), intent(inout) | ExtPtfm, | ||
type(map_data), intent(inout) | MAPp, | ||
type(feamooring_data), intent(inout) | FEAM, | ||
type(moordyn_data), intent(inout) | MD, | ||
type(orcaflex_data), intent(inout) | Orca, | ||
type(icefloe_data), intent(inout) | IceF, | ||
type(icedyn_data), intent(inout) | IceD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
logical, intent(in) | WriteThisStep | ||
) |
This routines advances the states of each module.
[in] | t_initial | initial simulation time (almost always 0) |
[in] | n_t_global | integer time step |
[in] | p_fast | Parameters for the glue code |
[in] | m_fast | Miscellaneous variables |
[in,out] | ed | ElastoDyn data |
[in,out] | bd | BeamDyn data |
[in,out] | srvd | ServoDyn data |
[in,out] | ad14 | AeroDyn v14 data |
[in,out] | ad | AeroDyn data |
[in,out] | ifw | InflowWind data |
[in,out] | opfm | OpenFOAM data |
[in,out] | hd | HydroDyn data |
[in,out] | sd | SubDyn data |
[in,out] | extptfm | ExtPtfm data |
[in,out] | mapp | MAP data |
[in,out] | feam | FEAMooring data |
[in,out] | md | Data for the MoorDyn module |
[in,out] | orca | OrcaFlex interface data |
[in,out] | icef | IceFloe data |
[in,out] | iced | All the IceDyn data used in time-step loop |
[in,out] | meshmapdata | Data for mapping between modules (added to help BD get better root motion inputs) |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
[in] | writethisstep | Will we print the WriteOutput values this step (for optimizations with SolveOption2)? |
subroutine fast_solver::fast_extrapinterpmods | ( | real(dbki), intent(in) | t_global_next, |
type(fast_parametertype), intent(in) | p_FAST, | ||
type(fast_miscvartype), intent(in) | m_FAST, | ||
type(elastodyn_data), intent(inout) | ED, | ||
type(beamdyn_data), intent(inout) | BD, | ||
type(servodyn_data), intent(inout) | SrvD, | ||
type(aerodyn14_data), intent(inout) | AD14, | ||
type(aerodyn_data), intent(inout) | AD, | ||
type(inflowwind_data), intent(inout) | IfW, | ||
type(hydrodyn_data), intent(inout) | HD, | ||
type(subdyn_data), intent(inout) | SD, | ||
type(extptfm_data), intent(inout) | ExtPtfm, | ||
type(map_data), intent(inout) | MAPp, | ||
type(feamooring_data), intent(inout) | FEAM, | ||
type(moordyn_data), intent(inout) | MD, | ||
type(orcaflex_data), intent(inout) | Orca, | ||
type(icefloe_data), intent(inout) | IceF, | ||
type(icedyn_data), intent(inout) | IceD, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine extrapolates inputs to modules to give predicted values at t+dt.
[in] | t_global_next | next global time step (t + dt), at which we're extrapolating inputs (and ED outputs) |
[in] | p_fast | Parameters for the glue code |
[in] | m_fast | Miscellaneous variables |
[in,out] | ed | ElastoDyn data |
[in,out] | bd | BeamDyn data |
[in,out] | srvd | ServoDyn data |
[in,out] | ad14 | AeroDyn14 data |
[in,out] | ad | AeroDyn data |
[in,out] | ifw | InflowWind data |
[in,out] | hd | HydroDyn data |
[in,out] | sd | SubDyn data |
[in,out] | extptfm | ExtPtfm data |
[in,out] | mapp | MAP data |
[in,out] | feam | FEAMooring data |
[in,out] | md | Data for the MoorDyn module |
[in,out] | orca | OrcaFlex interface data |
[in,out] | icef | IceFloe data |
[in,out] | iced | All the IceDyn data used in time-step loop |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::feam_inputsolve | ( | type(feam_inputtype), intent(inout) | u_FEAM, |
type(ed_outputtype), intent(in) | y_ED, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine sets the inputs required for FEAM.
[in,out] | u_feam | FEAM input |
[in] | y_ed | The outputs of the structural dynamics module |
[in,out] | meshmapdata | data for mapping meshes between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::fullopt1_inputoutputsolve | ( | real(dbki), intent(in) | this_time, |
type(fast_parametertype), intent(in) | p_FAST, | ||
logical, intent(in) | calcJacobian, | ||
type(ed_inputtype), intent(inout) | u_ED, | ||
type(ed_parametertype), intent(in) | p_ED, | ||
type(ed_continuousstatetype), intent(in) | x_ED, | ||
type(ed_discretestatetype), intent(in) | xd_ED, | ||
type(ed_constraintstatetype), intent(in) | z_ED, | ||
type(ed_otherstatetype), intent(in) | OtherSt_ED, | ||
type(ed_outputtype), intent(inout), target | y_ED, | ||
type(ed_miscvartype), intent(inout) | m_ED, | ||
type(sd_inputtype), intent(inout) | u_SD, | ||
type(sd_parametertype), intent(in) | p_SD, | ||
type(sd_continuousstatetype), intent(in) | x_SD, | ||
type(sd_discretestatetype), intent(in) | xd_SD, | ||
type(sd_constraintstatetype), intent(in) | z_SD, | ||
type(sd_otherstatetype), intent(in) | OtherSt_SD, | ||
type(sd_outputtype), intent(inout) | y_SD, | ||
type(sd_miscvartype), intent(inout) | m_SD, | ||
type(extptfm_inputtype), intent(inout) | u_ExtPtfm, | ||
type(extptfm_parametertype), intent(in) | p_ExtPtfm, | ||
type(extptfm_continuousstatetype), intent(in) | x_ExtPtfm, | ||
type(extptfm_discretestatetype), intent(in) | xd_ExtPtfm, | ||
type(extptfm_constraintstatetype), intent(in) | z_ExtPtfm, | ||
type(extptfm_otherstatetype), intent(in) | OtherSt_ExtPtfm, | ||
type(extptfm_outputtype), intent(inout) | y_ExtPtfm, | ||
type(extptfm_miscvartype), intent(inout) | m_ExtPtfm, | ||
type(hydrodyn_inputtype), intent(inout) | u_HD, | ||
type(hydrodyn_parametertype), intent(in) | p_HD, | ||
type(hydrodyn_continuousstatetype), intent(in) | x_HD, | ||
type(hydrodyn_discretestatetype), intent(in) | xd_HD, | ||
type(hydrodyn_constraintstatetype), intent(in) | z_HD, | ||
type(hydrodyn_otherstatetype), intent(inout) | OtherSt_HD, | ||
type(hydrodyn_outputtype), intent(inout) | y_HD, | ||
type(hydrodyn_miscvartype), intent(inout) | m_HD, | ||
type(bd_inputtype), dimension(:), intent(inout) | u_BD, | ||
type(bd_parametertype), dimension(:), intent(in) | p_BD, | ||
type(bd_continuousstatetype), dimension(:), intent(in) | x_BD, | ||
type(bd_discretestatetype), dimension(:), intent(in) | xd_BD, | ||
type(bd_constraintstatetype), dimension(:), intent(in) | z_BD, | ||
type(bd_otherstatetype), dimension(:), intent(in) | OtherSt_BD, | ||
type(bd_outputtype), dimension(:), intent(inout) | y_BD, | ||
type(bd_miscvartype), dimension(:), intent(inout) | m_BD, | ||
type(orca_inputtype), intent(inout) | u_Orca, | ||
type(orca_parametertype), intent(in) | p_Orca, | ||
type(orca_continuousstatetype), intent(in) | x_Orca, | ||
type(orca_discretestatetype), intent(in) | xd_Orca, | ||
type(orca_constraintstatetype), intent(in) | z_Orca, | ||
type(orca_otherstatetype), intent(in) | OtherSt_Orca, | ||
type(orca_outputtype), intent(inout) | y_Orca, | ||
type(orca_miscvartype), intent(inout) | m_Orca, | ||
type(map_inputtype), intent(inout) | u_MAP, | ||
type(map_outputtype), intent(in) | y_MAP, | ||
type(feam_inputtype), intent(inout) | u_FEAM, | ||
type(feam_outputtype), intent(in) | y_FEAM, | ||
type(md_inputtype), intent(inout) | u_MD, | ||
type(md_outputtype), intent(in) | y_MD, | ||
type(icefloe_inputtype), intent(inout) | u_IceF, | ||
type(icefloe_outputtype), intent(in) | y_IceF, | ||
type(iced_inputtype), dimension(:), intent(inout) | u_IceD, | ||
type(iced_outputtype), dimension(:), intent(in) | y_IceD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
logical, intent(in) | WriteThisStep | ||
) |
This routine performs the Input-Output solve for ED, SD, HD, BD, and/or the OrcaFlex Interface.
Note that this has been customized for the physics in the problems and is not a general solution.
[in] | this_time | The current simulation time (actual or time of prediction) |
[in] | p_fast | Glue-code simulation parameters |
[in] | calcjacobian | Should we calculate Jacobians this time? |
[in] | x_ed | Continuous states |
[in] | xd_ed | Discrete states |
[in] | z_ed | Constraint states |
[in] | otherst_ed | Other states |
[in] | p_ed | Parameters |
[in,out] | u_ed | System inputs |
[in,out] | y_ed | System outputs |
[in,out] | m_ed | misc/optimization variables |
[in] | x_bd | Continuous states |
[in] | xd_bd | Discrete states |
[in] | z_bd | Constraint states |
[in] | otherst_bd | Other/optimization states |
[in] | p_bd | Parameters |
[in,out] | u_bd | System inputs |
[in,out] | y_bd | System outputs |
[in,out] | m_bd | misc/optimization variables |
[in] | x_sd | Continuous states |
[in] | xd_sd | Discrete states |
[in] | z_sd | Constraint states |
[in] | otherst_sd | Other states |
[in] | p_sd | Parameters |
[in,out] | u_sd | System inputs |
[in,out] | y_sd | System outputs |
[in,out] | m_sd | misc/optimization variables |
[in] | x_extptfm | Continuous states |
[in] | xd_extptfm | Discrete states |
[in] | z_extptfm | Constraint states |
[in] | otherst_extptfm | Other states |
[in] | p_extptfm | Parameters |
[in,out] | u_extptfm | System inputs |
[in,out] | y_extptfm | System outputs |
[in,out] | m_extptfm | misc/optimization variables |
[in] | x_hd | Continuous states |
[in] | xd_hd | Discrete states |
[in] | z_hd | Constraint states |
[in,out] | otherst_hd | Other states |
[in] | p_hd | Parameters |
[in,out] | u_hd | System inputs |
[in,out] | y_hd | System outputs |
[in,out] | m_hd | misc/optimization variables |
[in] | x_orca | Continuous states |
[in] | xd_orca | Discrete states |
[in] | z_orca | Constraint states |
[in] | otherst_orca | Other states |
[in] | p_orca | Parameters |
[in,out] | u_orca | System inputs |
[in,out] | y_orca | System outputs |
[in,out] | m_orca | misc/optimization variables |
[in] | y_map | MAP outputs |
[in,out] | u_map | MAP inputs (INOUT just because I don't want to use another tempoarary mesh and we'll overwrite this later) |
[in] | y_feam | FEAM outputs |
[in,out] | u_feam | FEAM inputs (INOUT just because I don't want to use another tempoarary mesh and we'll overwrite this later) |
[in] | y_md | MoorDyn outputs |
[in,out] | u_md | MoorDyn inputs (INOUT just because I don't want to use another tempoarary mesh and we'll overwrite this later) |
[in] | y_icef | IceFloe outputs |
[in,out] | u_icef | IceFloe inputs (INOUT just because I don't want to use another tempoarary mesh and we'll overwrite this later) |
[in] | y_iced | IceDyn outputs |
[in,out] | u_iced | IceDyn inputs (INOUT just because I don't want to use another tempoarary mesh and we'll overwrite this later) |
[in,out] | meshmapdata | data for mapping meshes between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
[in] | writethisstep | Will we print the WriteOutput values this step? |
real(reki) function fast_solver::getperturb | ( | real(reki), intent(in) | x | ) |
function to return the size of perturbation in calculating jacobian with finite differences.
Currently hard-coded to return 1.
[in] | x | value that we want to perturb |
subroutine fast_solver::iced_inputsolve | ( | type(iced_inputtype), intent(inout) | u_IceD, |
type(sd_outputtype), intent(in) | y_SD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(in) | legNum, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine sets the inputs required for IceFloe.
[in,out] | u_iced | IceDyn input |
[in] | y_sd | SubDyn outputs |
[in,out] | meshmapdata | data for mapping meshes between modules |
[in] | legnum | which instance of IceDyn we're using |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::icefloe_inputsolve | ( | type(icefloe_inputtype), intent(inout) | u_IceF, |
type(sd_outputtype), intent(in) | y_SD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine sets the inputs required for IceFloe.
[in,out] | u_icef | IceFloe input |
[in] | y_sd | SubDyn outputs |
[in,out] | meshmapdata | data for mapping meshes between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::ifw_inputsolve | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(fast_miscvartype), intent(in) | m_FAST, | ||
type(inflowwind_inputtype), dimension(:), intent(inout) | u_IfW, | ||
type(inflowwind_parametertype), intent(in) | p_IfW, | ||
type(ad14_inputtype), intent(in) | u_AD14, | ||
type(ad_inputtype), intent(in) | u_AD, | ||
type(ad_otherstatetype), intent(in) | OtherSt_AD, | ||
type(ed_outputtype), intent(in) | y_ED, | ||
integer(intki) | ErrStat, | ||
character(*) | ErrMsg | ||
) |
This routine determines the points in space where InflowWind needs to compute wind speeds.
[in,out] | u_ifw | The inputs to InflowWind |
[in] | p_ifw | The parameters to InflowWind |
[in] | u_ad14 | The input meshes (already calculated) from AeroDyn14 |
[in] | u_ad | The input meshes (already calculated) from AeroDyn |
[in] | otherst_ad | The wake points from AeroDyn are in here (Free Vortex Wake) |
[in] | y_ed | The outputs of the structural dynamics module (for IfW Lidar) |
[in] | p_fast | FAST parameter data |
[in] | m_fast | misc FAST data, including inputs from external codes like Simulink |
errstat | Error status of the operation | |
errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::ifw_setexternalinputs | ( | type(inflowwind_parametertype), intent(in) | p_IfW, |
type(fast_miscvartype), intent(in) | m_FAST, | ||
type(ed_outputtype), intent(in) | y_ED, | ||
type(inflowwind_inputtype), intent(inout) | u_IfW | ||
) |
This routine sets the inputs required for InflowWind from an external source (Simulink)
[in] | p_ifw | InflowWind parameters |
[in] | m_fast | Glue-code misc variables (including inputs from external sources like Simulink) |
[in] | y_ed | The outputs of the structural dynamics module |
[in,out] | u_ifw | InflowWind Inputs at t |
subroutine fast_solver::init_fullopt1_jacobian | ( | type(fast_parametertype), intent(inout) | p_FAST, |
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
type(meshtype), intent(in) | ED_PlatformPtMesh, | ||
type(meshtype), intent(in) | SD_TPMesh, | ||
type(meshtype), intent(in) | SD_LMesh, | ||
type(meshtype), intent(in) | HD_M_LumpedMesh, | ||
type(meshtype), intent(in) | HD_M_DistribMesh, | ||
type(meshtype), intent(in) | HD_WAMIT_Mesh, | ||
type(meshtype), intent(in) | ED_HubPtLoad, | ||
type(bd_inputtype), dimension(:), intent(in) | u_BD, | ||
type(meshtype), intent(in) | Orca_PtfmMesh, | ||
type(meshtype), intent(in) | ExtPtfm_PtfmMesh, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine initializes the array that maps rows/columns of the Jacobian to specific mesh fields.
Do not change the order of this packing without changing subroutine Create_FullOpt1_UVector()!
[in,out] | p_fast | FAST parameters |
[in,out] | meshmapdata | data that maps meshes together |
[in] | ed_platformptmesh | ElastoDyn's PlatformPtMesh |
[in] | ed_hubptload | ElastoDyn's HubPtLoad mesh |
[in] | sd_tpmesh | SubDyn's TP (transition piece) mesh |
[in] | sd_lmesh | SubDyn's LMesh |
[in] | hd_m_lumpedmesh | HydroDyn's Morison Lumped Mesh |
[in] | hd_m_distribmesh | HydroDyn's Morison Distributed Mesh |
[in] | hd_wamit_mesh | HydroDyn's WAMIT mesh |
[in] | u_bd | inputs for each instance of the BeamDyn module (for the RootMotion meshes) |
[in] | orca_ptfmmesh | OrcaFlex interface PtfmMesh |
[in] | extptfm_ptfmmesh | ExtPtfm_MCKF interface PtfmMesh |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::initmodulemappings | ( | type(fast_parametertype), intent(inout) | p_FAST, |
type(elastodyn_data), intent(inout), target | ED, | ||
type(beamdyn_data), intent(inout) | BD, | ||
type(aerodyn14_data), intent(inout) | AD14, | ||
type(aerodyn_data), intent(inout) | AD, | ||
type(hydrodyn_data), intent(inout) | HD, | ||
type(subdyn_data), intent(inout) | SD, | ||
type(extptfm_data), intent(inout) | ExtPtfm, | ||
type(servodyn_data), intent(inout) | SrvD, | ||
type(map_data), intent(inout) | MAPp, | ||
type(feamooring_data), intent(inout) | FEAM, | ||
type(moordyn_data), intent(inout) | MD, | ||
type(orcaflex_data), intent(inout) | Orca, | ||
type(icefloe_data), intent(inout) | IceF, | ||
type(icedyn_data), intent(inout) | IceD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine initializes all of the mapping data structures needed between the various modules.
[in,out] | p_fast | Parameters for the glue code |
[in,out] | ed | ElastoDyn data |
[in,out] | bd | BeamDyn data |
[in,out] | srvd | ServoDyn data |
[in,out] | ad | AeroDyn data |
[in,out] | ad14 | AeroDyn14 data |
[in,out] | hd | HydroDyn data |
[in,out] | sd | SubDyn data |
[in,out] | extptfm | ExtPtfm data |
[in,out] | mapp | MAP data |
[in,out] | feam | FEAMooring data |
[in,out] | md | MoorDyn data |
[in,out] | orca | OrcaFlex interface data |
[in,out] | icef | IceFloe data |
[in,out] | iced | All the IceDyn data used in time-step loop |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::map_inputsolve | ( | type(map_inputtype), intent(inout) | u_MAP, |
type(ed_outputtype), intent(in) | y_ED, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine sets the inputs required for MAP.
[in,out] | u_map | MAP input |
[in] | y_ed | The outputs of the structural dynamics module |
[in,out] | meshmapdata | data for mapping meshes between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::md_inputsolve | ( | type(md_inputtype), intent(inout) | u_MD, |
type(ed_outputtype), intent(in) | y_ED, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine sets the inputs required for MoorDyn.
[in,out] | u_md | MoorDyn input |
[in] | y_ed | The outputs of the structural dynamics module |
[in,out] | meshmapdata | data for mapping meshes between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::perturb_u_fullopt1 | ( | type(fast_parametertype), intent(in) | p_FAST, |
integer( intki ), dimension(:,:), intent(in) | Jac_u_indx, | ||
integer( intki ), intent(in) | n, | ||
real( reki ), dimension(:), intent(inout) | u_perturb, | ||
type(ed_inputtype), intent(inout), optional | u_ED_perturb, | ||
type(sd_inputtype), intent(inout), optional | u_SD_perturb, | ||
type(hydrodyn_inputtype), intent(inout), optional | u_HD_perturb, | ||
type(bd_inputtype), intent(inout), optional | u_BD_perturb, | ||
type(orca_inputtype), intent(inout), optional | u_Orca_perturb, | ||
type(extptfm_inputtype), intent(inout), optional | u_ExtPtfm_perturb, | ||
real( reki ), intent(out) | perturb | ||
) |
This routine perturbs the nth element of the u array (and mesh/field it corresponds to)
[in] | p_fast | Glue-code simulation parameters |
[in] | jac_u_indx | Index to map Jacobian u-vector into mesh fields |
[in] | n | number of array element to use |
[in,out] | u_perturb | array to be perturbed |
[in,out] | u_ed_perturb | ED System inputs (needed only when 1 <= n <= NumEDNodes=NumEDNodes) |
[in,out] | u_sd_perturb | SD System inputs (needed only when NumEDNodes +1 <= n <= NumEDNodes+NumSDNodes) [if SD is used] |
[in,out] | u_hd_perturb | HD System inputs (needed only when NumEDNodes+NumSDNodes +1 <= n <= NumEDNodes+NumSDNodes+NumHDNodes) [if HD is used and SD is used. if SD not used, |
[in,out] | u_bd_perturb | BD System inputs (needed only when NumEDNodes+NumSDNodes+NumHDNodes+1 <= n <= inf) [if BD is used] |
[in,out] | u_orca_perturb | Orca System inputs (needed only when NumEDNodes+NumSDNodes+NumHDNodes+NumBDNodes+1 <= n <= inf) [if Orca is used] |
[in,out] | u_extptfm_perturb | ExtPtfm System inputs (needed only when NumEDNodes+NumSDNodes+NumHDNodes+NumBDNodes+NumOcraNodes+1 <= n <= inf) [if ExtPtfm is used] |
[out] | perturb | amount that u_perturb(n) was perturbed |
subroutine fast_solver::resetremapflags | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(elastodyn_data), intent(inout) | ED, | ||
type(beamdyn_data), intent(inout) | BD, | ||
type(aerodyn14_data), intent(inout) | AD14, | ||
type(aerodyn_data), intent(inout) | AD, | ||
type(hydrodyn_data), intent(inout) | HD, | ||
type(subdyn_data), intent(inout) | SD, | ||
type(extptfm_data), intent(inout) | ExtPtfm, | ||
type(servodyn_data), intent(inout) | SrvD, | ||
type(map_data), intent(inout) | MAPp, | ||
type(feamooring_data), intent(inout) | FEAM, | ||
type(moordyn_data), intent(inout) | MD, | ||
type(orcaflex_data), intent(inout) | Orca, | ||
type(icefloe_data), intent(inout) | IceF, | ||
type(icedyn_data), intent(inout) | IceD | ||
) |
This routine resets the remap flags on all of the meshes.
[in] | p_fast | Parameters for the glue code |
[in,out] | ed | ElastoDyn data |
[in,out] | bd | BeamDyn data |
[in,out] | srvd | ServoDyn data |
[in,out] | ad | AeroDyn data |
[in,out] | ad14 | AeroDyn14 data |
[in,out] | hd | HydroDyn data |
[in,out] | sd | SubDyn data |
[in,out] | extptfm | ExtPtfm data |
[in,out] | mapp | MAP data |
[in,out] | feam | FEAMooring data |
[in,out] | md | MoorDyn data |
[in,out] | orca | OrcaFlex interface data |
[in,out] | icef | IceFloe data |
[in,out] | iced | All the IceDyn data used in time-step loop |
subroutine fast_solver::solveoption1 | ( | real(dbki), intent(in) | this_time, |
integer(intki), intent(in) | this_state, | ||
logical, intent(in) | calcJacobian, | ||
type(fast_parametertype), intent(in) | p_FAST, | ||
type(elastodyn_data), intent(inout) | ED, | ||
type(beamdyn_data), intent(inout) | BD, | ||
type(hydrodyn_data), intent(inout) | HD, | ||
type(subdyn_data), intent(inout) | SD, | ||
type(extptfm_data), intent(inout) | ExtPtfm, | ||
type(map_data), intent(inout) | MAPp, | ||
type(feamooring_data), intent(inout) | FEAM, | ||
type(moordyn_data), intent(inout) | MD, | ||
type(orcaflex_data), intent(inout) | Orca, | ||
type(icefloe_data), intent(inout) | IceF, | ||
type(icedyn_data), intent(inout) | IceD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
logical, intent(in) | WriteThisStep | ||
) |
This routine implements the "option 1" solve for all inputs with direct links to HD, SD, ExtPtfm, MAP, OrcaFlex interface, and the ED platform reference point.
Also in solve option 1 are the BD-ED blade root coupling.
[in] | this_time | The current simulation time (actual or time of prediction) |
[in] | this_state | Index into the state array (current or predicted states) |
[in] | calcjacobian | Should we calculate Jacobians in Option 1? |
[in] | p_fast | Parameters for the glue code |
[in,out] | ed | ElastoDyn data |
[in,out] | bd | BeamDyn data |
[in,out] | hd | HydroDyn data |
[in,out] | sd | SubDyn data |
[in,out] | extptfm | ExtPtfm data |
[in,out] | mapp | MAP data |
[in,out] | feam | FEAMooring data |
[in,out] | md | MoorDyn data |
[in,out] | orca | OrcaFlex interface data |
[in,out] | icef | IceFloe data |
[in,out] | iced | All the IceDyn data used in time-step loop |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
[in] | writethisstep | Will we print the WriteOutput values this step? |
Option 1: solve for consistent inputs and outputs, which is required when Y has direct feedthrough in modules coupled together
subroutine fast_solver::solveoption2 | ( | real(dbki), intent(in) | this_time, |
integer(intki), intent(in) | this_state, | ||
type(fast_parametertype), intent(in) | p_FAST, | ||
type(fast_miscvartype), intent(in) | m_FAST, | ||
type(elastodyn_data), intent(inout) | ED, | ||
type(beamdyn_data), intent(inout) | BD, | ||
type(aerodyn14_data), intent(inout) | AD14, | ||
type(aerodyn_data), intent(inout) | AD, | ||
type(servodyn_data), intent(inout) | SrvD, | ||
type(inflowwind_data), intent(inout) | IfW, | ||
type(openfoam_data), intent(inout) | OpFM, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
logical, intent(in) | firstCall, | ||
logical, intent(in) | WriteThisStep | ||
) |
This routine implements the "option 2" solve for all inputs without direct links to HD, SD, MAP, or the ED platform reference point.
[in] | firstcall | flag to determine how to call ServoDyn (a hack) |
[in] | this_time | The current simulation time (actual or time of prediction) |
[in] | this_state | Index into the state array (current or predicted states) |
[in] | p_fast | Parameters for the glue code |
[in] | m_fast | Misc variables for the glue code (including external inputs) |
[in,out] | ed | ElastoDyn data |
[in,out] | bd | BeamDyn data |
[in,out] | srvd | ServoDyn data |
[in,out] | ad14 | AeroDyn14 data |
[in,out] | ad | AeroDyn data |
[in,out] | ifw | InflowWind data |
[in,out] | opfm | OpenFOAM data |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
[in] | writethisstep | Will we print the WriteOutput values this step? |
++ Option 2: Solve for inputs based only on the current outputs. This is much faster than option 1 when the coupled modules do not have direct feedthrough.
subroutine fast_solver::solveoption2a_inp2bd | ( | real(dbki), intent(in) | this_time, |
integer(intki), intent(in) | this_state, | ||
type(fast_parametertype), intent(in) | p_FAST, | ||
type(fast_miscvartype), intent(in) | m_FAST, | ||
type(elastodyn_data), intent(inout) | ED, | ||
type(beamdyn_data), intent(inout) | BD, | ||
type(aerodyn14_data), intent(inout) | AD14, | ||
type(aerodyn_data), intent(inout) | AD, | ||
type(servodyn_data), intent(inout) | SrvD, | ||
type(inflowwind_data), intent(inout) | IfW, | ||
type(openfoam_data), intent(inout) | OpFM, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
logical, intent(in) | WriteThisStep | ||
) |
This routine implements the first part of the "option 2" solve for inputs that apply to BeamDyn and AeroDyn.
[in] | this_time | The current simulation time (actual or time of prediction) |
[in] | this_state | Index into the state array (current or predicted states) |
[in] | p_fast | Parameters for the glue code |
[in] | m_fast | Misc variables for the glue code (including external inputs) |
[in,out] | ed | ElastoDyn data |
[in,out] | bd | BeamDyn data |
[in,out] | srvd | ServoDyn data |
[in,out] | ad14 | AeroDyn14 data |
[in,out] | ad | AeroDyn data |
[in,out] | ifw | InflowWind data |
[in,out] | opfm | OpenFOAM data |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
[in] | writethisstep | Will we print the WriteOutput values this step? |
++ Option 2: Solve for inputs based only on the current outputs. This is much faster than option 1 when the coupled modules do not have direct feedthrough.
subroutine fast_solver::solveoption2b_inp2ifw | ( | real(dbki), intent(in) | this_time, |
integer(intki), intent(in) | this_state, | ||
type(fast_parametertype), intent(in) | p_FAST, | ||
type(fast_miscvartype), intent(in) | m_FAST, | ||
type(elastodyn_data), intent(inout) | ED, | ||
type(beamdyn_data), intent(inout) | BD, | ||
type(aerodyn14_data), intent(inout) | AD14, | ||
type(aerodyn_data), intent(inout) | AD, | ||
type(servodyn_data), intent(inout) | SrvD, | ||
type(inflowwind_data), intent(inout) | IfW, | ||
type(openfoam_data), intent(inout) | OpFM, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
logical, intent(in) | WriteThisStep | ||
) |
This routine implements the first part of the "option 2" solve for inputs that apply to AeroDyn & InflowWind.
[in] | this_time | The current simulation time (actual or time of prediction) |
[in] | this_state | Index into the state array (current or predicted states) |
[in] | p_fast | Parameters for the glue code |
[in] | m_fast | Misc variables for the glue code (including external inputs) |
[in,out] | ed | ElastoDyn data |
[in,out] | bd | BeamDyn data |
[in,out] | srvd | ServoDyn data |
[in,out] | ad14 | AeroDyn14 data |
[in,out] | ad | AeroDyn data |
[in,out] | ifw | InflowWind data |
[in,out] | opfm | OpenFOAM data |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
[in] | writethisstep | Will we print the WriteOutput values this step? |
++ Option 2: Solve for inputs based only on the current outputs. This is much faster than option 1 when the coupled modules do not have direct feedthrough.
subroutine fast_solver::solveoption2c_inp2ad_srvd | ( | real(dbki), intent(in) | this_time, |
integer(intki), intent(in) | this_state, | ||
type(fast_parametertype), intent(in) | p_FAST, | ||
type(fast_miscvartype), intent(in) | m_FAST, | ||
type(elastodyn_data), intent(inout) | ED, | ||
type(beamdyn_data), intent(inout) | BD, | ||
type(aerodyn14_data), intent(inout) | AD14, | ||
type(aerodyn_data), intent(inout) | AD, | ||
type(servodyn_data), intent(inout) | SrvD, | ||
type(inflowwind_data), intent(inout) | IfW, | ||
type(openfoam_data), intent(inout) | OpFM, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
logical, intent(in) | WriteThisStep | ||
) |
This routine implements the first part of the "option 2" solve for inputs that apply to AeroDyn and ServoDyn.
[in] | this_time | The current simulation time (actual or time of prediction) |
[in] | this_state | Index into the state array (current or predicted states) |
[in] | p_fast | Parameters for the glue code |
[in] | m_fast | Misc variables for the glue code (including external inputs) |
[in,out] | ed | ElastoDyn data |
[in,out] | bd | BeamDyn data |
[in,out] | srvd | ServoDyn data |
[in,out] | ad14 | AeroDyn14 data |
[in,out] | ad | AeroDyn data |
[in,out] | ifw | InflowWind data |
[in,out] | opfm | OpenFOAM data |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
[in] | writethisstep | Will we print the WriteOutput values this step? |
++ Option 2: Solve for inputs based only on the current outputs. This is much faster than option 1 when the coupled modules do not have direct feedthrough.
subroutine fast_solver::srvd_inputsolve | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(fast_miscvartype), intent(in) | m_FAST, | ||
type(srvd_inputtype), intent(inout) | u_SrvD, | ||
type(ed_outputtype), intent(in) | y_ED, | ||
type(inflowwind_outputtype), intent(in) | y_IfW, | ||
type(opfm_outputtype), intent(in) | y_OpFM, | ||
type(bd_outputtype), dimension(:), intent(in) | y_BD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine sets the inputs required for ServoDyn.
[in] | p_fast | Glue-code simulation parameters |
[in] | m_fast | Glue-code misc variables (including inputs from external sources like Simulink) |
[in,out] | u_srvd | ServoDyn Inputs at t |
[in] | y_ed | ElastoDyn outputs |
[in] | y_ifw | InflowWind outputs |
[in] | y_opfm | OpenFOAM outputs |
[in] | y_bd | BD Outputs |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status |
[out] | errmsg | Error message |
subroutine fast_solver::srvd_setexternalinputs | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(fast_miscvartype), intent(in) | m_FAST, | ||
type(srvd_inputtype), intent(inout) | u_SrvD | ||
) |
This routine sets the inputs required for ServoDyn from an external source (Simulink)
[in] | p_fast | Glue-code simulation parameters |
[in] | m_fast | Glue-code misc variables (including inputs from external sources like Simulink) |
[in,out] | u_srvd | ServoDyn Inputs at t |
subroutine fast_solver::transfer_ed_to_bd | ( | type(ed_outputtype), intent(in) | y_ED, |
type(bd_inputtype), dimension(:), intent(inout) | u_BD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine sets the inputs required for IceFloe.
[in,out] | u_bd | BeamDyn inputs |
[in] | y_ed | ElastoDyn outputs |
[in,out] | meshmapdata | data for mapping meshes between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::transfer_ed_to_bd_tmp | ( | type(ed_outputtype), intent(in) | y_ED, |
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine sets the inputs required for IceFloe.
[in] | y_ed | ElastoDyn outputs |
[in,out] | meshmapdata | data for mapping meshes between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::transfer_ed_to_hd_sd_bd_mooring | ( | type(fast_parametertype), intent(in) | p_FAST, |
type(ed_outputtype), intent(in) | y_ED, | ||
type(hydrodyn_inputtype), intent(inout) | u_HD, | ||
type(sd_inputtype), intent(inout) | u_SD, | ||
type(extptfm_inputtype), intent(inout) | u_ExtPtfm, | ||
type(map_inputtype), intent(inout) | u_MAP, | ||
type(feam_inputtype), intent(inout) | u_FEAM, | ||
type(md_inputtype), intent(inout) | u_MD, | ||
type(orca_inputtype), intent(inout) | u_Orca, | ||
type(bd_inputtype), dimension(:), intent(inout) | u_BD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine transfers the ED outputs into inputs required for HD, SD, ExtPtfm, BD, MAP, and/or FEAM.
[in] | p_fast | Glue-code simulation parameters |
[in] | y_ed | The outputs of the structural dynamics module |
[in,out] | u_hd | HydroDyn input |
[in,out] | u_sd | SubDyn input |
[in,out] | u_extptfm | ExtPtfm_MCKF input |
[in,out] | u_map | MAP input |
[in,out] | u_feam | FEAM input |
[in,out] | u_md | MoorDyn input |
[in,out] | u_orca | OrcaFlex input |
[in,out] | u_bd | BeamDyn inputs |
[in,out] | meshmapdata | data for mapping meshes between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::transfer_hd_to_sd | ( | type(meshtype), intent(inout) | u_mapped, |
type(meshtype), intent(inout) | u_SD_LMesh, | ||
type(meshtype), intent(in) | u_mapped_positions, | ||
type(hydrodyn_outputtype), intent(in) | y_HD, | ||
type(meshtype), intent(in) | u_HD_M_LumpedMesh, | ||
type(meshtype), intent(in) | u_HD_M_DistribMesh, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine transfers the HD outputs into inputs required for ED.
Note that this adds to the values already in u_SD_LMesh (so initialize it before calling this routine).
[in,out] | u_mapped | temporary copy of SD mesh (an argument to avoid another temporary mesh copy) |
[in,out] | u_sd_lmesh | SD Inputs on LMesh at t (separate so we can call from FullOpt1_InputOutputSolve with temp meshes) |
[in] | u_mapped_positions | Mesh sibling of u_mapped, with displaced positions |
[in] | y_hd | HydroDyn outputs |
[in] | u_hd_m_lumpedmesh | HydroDyn input mesh (separate so we can call from FullOpt1_InputOutputSolve with temp meshes) |
[in] | u_hd_m_distribmesh | HydroDyn input mesh (separate so we can call from FullOpt1_InputOutputSolve with temp meshes) |
[in,out] | meshmapdata | Data for mapping between modules |
[out] | errstat | Error status |
[out] | errmsg | Error message |
subroutine fast_solver::transfer_platformmotion_to_hd | ( | type(meshtype), intent(in) | PlatformMotion, |
type(hydrodyn_inputtype), intent(inout) | u_HD, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine transfers the platform motion output of the structural module (ED) into inputs required for HD.
[in] | platformmotion | The platform motion outputs of the structural dynamics module |
[in,out] | u_hd | HydroDyn input |
[in,out] | meshmapdata | data for mapping meshes between modules |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine fast_solver::transfer_sd_to_hd | ( | type(sd_outputtype), intent(in) | y_SD, |
type(meshtype), intent(inout) | u_HD_M_LumpedMesh, | ||
type(meshtype), intent(inout) | u_HD_M_DistribMesh, | ||
type(fast_modulemaptype), intent(inout) | MeshMapData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine transfers the SD outputs into inputs required for HD.
[in] | y_sd | The outputs of the structural dynamics module |
[in,out] | u_hd_m_lumpedmesh | HydroDyn input mesh (separated here so that we can use temp meshes in ED_SD_HD_InputSolve) |
[in,out] | u_hd_m_distribmesh | HydroDyn input mesh (separated here so that we can use temp meshes in ED_SD_HD_InputSolve) |
[in,out] | meshmapdata | data for mapping meshes |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |