From abbb0741621d1871067e6ebee17a83493b55d0ab Mon Sep 17 00:00:00 2001 From: rebeccamccabe Date: Sat, 30 Nov 2024 15:58:25 -0500 Subject: [PATCH 1/7] spar buckling runs when set to compact --- mdocean/inputs/parameters.m | 1 + mdocean/simulation/modules/structures.m | 64 +++++++++++++++++++++---- mdocean/simulation/simulation.m | 2 +- 3 files changed, 58 insertions(+), 9 deletions(-) diff --git a/mdocean/inputs/parameters.m b/mdocean/inputs/parameters.m index fe7b0028..95eee7cb 100644 --- a/mdocean/inputs/parameters.m +++ b/mdocean/inputs/parameters.m @@ -61,6 +61,7 @@ ...% RM3 CBS sheet 1.4 average of cells F21, F34, F46, F56 ...% https://www.concretenetwork.com/concrete-prices.html ...% https://agmetalminer.com/metal-prices/ + table("nu","\nu",0.3,"structures",true,"Poisson's ratio (-)"); ... ...% Thicknesses and structures table("t_f_t","t_{f,t}",{0.50 * in2m},"structures",true,"float top thickness (m)"); diff --git a/mdocean/simulation/modules/structures.m b/mdocean/simulation/modules/structures.m index 47ef217a..e1d0d2df 100644 --- a/mdocean/simulation/modules/structures.m +++ b/mdocean/simulation/modules/structures.m @@ -1,6 +1,6 @@ function [FOS1Y,FOS2Y,FOS3Y,FOS_buckling] = structures(... F_heave, F_surge, M, h_s, T_s, rho_w, g, ... - sigma_y, A_c, A_lat_sub, r_over_t, I, E) + sigma_y, A_c, A_lat_sub, D_s, t_s_r, I, E, nu) %% Stress calculations depth = [0 T_s T_s]; % max depth @@ -9,7 +9,7 @@ sigma_surge = F_surge ./ A_lat_sub; sigma_rr = P_hydrostatic + sigma_surge; % radial compression - sigma_tt = P_hydrostatic .* r_over_t; % hoop stress + sigma_tt = 0;%P_hydrostatic .* r_over_t; % hoop stress sigma_zz = F_heave ./ A_c; % axial compression sigma_rt = sigma_surge; % shear sigma_tz = [0 0 0]; @@ -27,19 +27,67 @@ sigma_vm = von_mises(sigma_rr, sigma_tt, sigma_zz, sigma_rt, sigma_tz, sigma_zr); %% Buckling calculation - K = 2; % fixed-free - top is fixed by float angular stiffness, bottom is free - L = h_s; - F_buckling = pi^2 * E(M) * I(2) / (K*L)^2; - + FOS_spar = spar_combined_buckling(F_heave, E(M), I(2), h_s, D_s, A_c(2), t_s_r, ... + P_hydrostatic(2), sigma_y(M), nu); + %% Factor of Safety (FOS) Calculations FOS_yield = sigma_y(M) ./ sigma_vm; FOS1Y = FOS_yield(1); - FOS2Y = FOS_yield(2); + FOS2Y = FOS_spar; FOS3Y = FOS_yield(3); - FOS_buckling = F_buckling ./ F_heave; + + FOS_buckling = 3; end +function [FOS_spar, FOS_spar_local] = spar_combined_buckling(F, E, I, L, D, A, t, q, sigma_0, nu) + % euler buckling + K = 2; % fixed-free - top is fixed by float angular stiffness, bottom is free + F_buckling = pi^2 * E * I / (K*L)^2; + sigma_EA = F_buckling / A; + + % hoop stress + sigma_theta = q*D/(2*t); + + % local buckling of plate element: 2/9.7 + k_s = 1.33; % uniform compression, fixed-free, from table 3 on page 30 + s = D; % not sure if this is correct + P_r = 0.6; % steel proportional linear elastic limit + sigma_Ex = k_s * pi^2 * E / (12 * (1-nu^2)) * (t/s)^2; + if sigma_Ex <= P_r * sigma_0 + sigma_Cx = sigma_Ex; + else + sigma_Cx = sigma_0 * (1 - P_r * (1-P_r) * sigma_0/sigma_Ex); + end + + % stress at failure: 2/7.3 + compact = true; % assumption, eventually this should be a calculation + if compact + sigma_F = sigma_0; % yield + else + sigma_F = sigma_Cx; % local buckling of plate element + end + + % whether the failure mode is pure euler buckling or combined loading + sigma_EA_thresh = P_r * sigma_F * (1 - sigma_theta / sigma_F); + pure_buckling = sigma_EA <= sigma_EA_thresh; + if pure_buckling + sigma_C_A_theta = sigma_EA; + else + zeta = 1 - P_r*(1 - P_r)* sigma_F/sigma_EA - sigma_theta/sigma_F; + omega = 0.5 * sigma_theta/sigma_F * (1 - .5 * sigma_theta/sigma_F); + Lambda = 1/2 * (zeta + sqrt(zeta^2 + 4*omega)); + sigma_C_A_theta = sigma_F * Lambda; + end + + sigma_ac = F ./ A + q; + FOS_spar = sigma_C_A_theta / sigma_ac; + + % local buckling of tube: final part of 2/7.3, which combines 2/9.1 + % and 2/9.5 + FOS_spar_local = 0; +end + function s_vm = von_mises(s_11, s_22, s_33, s_12, s_23, s_31) principal_term = 1/2 * ( (s_11 - s_22).^2 + (s_22 - s_33).^2 + (s_33 - s_11).^2 ); diff --git a/mdocean/simulation/simulation.m b/mdocean/simulation/simulation.m index 6f954513..522a5dee 100644 --- a/mdocean/simulation/simulation.m +++ b/mdocean/simulation/simulation.m @@ -50,7 +50,7 @@ [FOS1Y,FOS2Y,FOS3Y,FOS_buckling] = structures(... F_heave_max, F_surge_max,... in.M, in.h_s, in.T_s, in.rho_w, in.g, in.sigma_y, A_c, ... - A_lat_sub, r_over_t, I, in.E); + A_lat_sub, in.D_s, in.t_s_r, I, in.E, in.nu); LCOE = econ(m_m, in.M, in.cost_m, in.N_WEC, P_avg_elec, in.FCR, in.cost_perN, in.F_max, in.eff_array); From 3cd7c365cc308613894511c000263fa21930d200 Mon Sep 17 00:00:00 2001 From: rebeccamccabe Date: Sun, 1 Dec 2024 00:32:15 -0500 Subject: [PATCH 2/7] more descriptive FOS names, and use FOS=11 not 3 for report validation --- mdocean/inputs/validation/validation_inputs.m | 4 +++- mdocean/simulation/simulation.m | 15 ++++++++------- 2 files changed, 11 insertions(+), 8 deletions(-) diff --git a/mdocean/inputs/validation/validation_inputs.m b/mdocean/inputs/validation/validation_inputs.m index 59e874d7..3e47ab56 100644 --- a/mdocean/inputs/validation/validation_inputs.m +++ b/mdocean/inputs/validation/validation_inputs.m @@ -17,7 +17,9 @@ 'power_avg', 85.9e3, ... % S14 in CBS Performance & Economics 'power_max', 286e3, ... % S15 in CBS Performance & Economics 'force_heave', 8500e3,... % p156 RM3 report - 'FOS_b', 3,... % p158 RM3 report + 'FOS_spar', 11.1,... % p158 RM3 report, but modified + ... % assuming that RM3 report used r + ... % instead of D when calculating I 'c_v', nominal_c_v()); % calculated from data in CBS Performance & Economics RM3_wecsim = struct(... diff --git a/mdocean/simulation/simulation.m b/mdocean/simulation/simulation.m index 522a5dee..c5181d69 100644 --- a/mdocean/simulation/simulation.m +++ b/mdocean/simulation/simulation.m @@ -47,7 +47,8 @@ P_var, P_avg_elec, P_matrix_elec, ... X_constraints] = dynamics(in, m_f_tot, m_s_tot, V_d, T); -[FOS1Y,FOS2Y,FOS3Y,FOS_buckling] = structures(... +[FOS_float,FOS_spar,FOS_damping_plate,... + FOS_spar_local_buckling] = structures(... F_heave_max, F_surge_max,... in.M, in.h_s, in.T_s, in.rho_w, in.g, in.sigma_y, A_c, ... A_lat_sub, in.D_s, in.t_s_r, I, in.E, in.nu); @@ -62,10 +63,10 @@ g(3) = V_s_pct; % prevent spar too heavy g(4) = 1 - V_s_pct; % prevent spar too light g(5) = GM; % pitch stability of float-spar system -g(6) = FOS1Y / p.FOS_min - 1; % float survives max force -g(7) = FOS2Y / p.FOS_min - 1; % spar survives max force -g(8) = FOS3Y / p.FOS_min - 1; % damping plate survives max force -g(9) = FOS_buckling / p.FOS_min - 1; % spar survives max force in buckling +g(6) = FOS_float / p.FOS_min - 1; % float survives max force +g(7) = FOS_spar / p.FOS_min - 1; % spar survives max force +g(8) = FOS_damping_plate / p.FOS_min - 1; % damping plate survives max force +g(9) = FOS_spar_local_buckling / p.FOS_min - 1; % spar survives max force in buckling g(10) = P_avg_elec; % positive power %1 + min(Kp_over_Ks,[],'all'); % spar heave stability (positive effective stiffness) g(11) = p.LCOE_max/LCOE - 1; % prevent more expensive than threshold @@ -85,7 +86,7 @@ end if nargout > 4 % if returning extra struct output for validation - [~, ~, ~, ~, ~, ~, ~, ~, ~, ~, ~, ~, ... + [~, ~, ~, ~, ~, ~, ~, ~, ~, ~, ~, ... mass, CB_f,CG_f] = geometry(in.D_s, in.D_f, in.D_f_in, in.D_f_b, ... in.T_f_1, in.T_f_2, in.h_f, in.h_s, ... in.h_fs_clear, in.D_f_tu, in.t_f_t, ... @@ -106,7 +107,7 @@ val.power_max = max(P_matrix_elec,[],'all'); val.force_heave = F_heave_max; val.force_ptrain = F_ptrain_max; - val.FOS_b = FOS_buckling; + val.FOS_spar = FOS_spar; val.c_v = P_var; val.B_p = B_p; val.X_u = X_u; From a756a40e777fbaf41616efa983927cf5ebb76062 Mon Sep 17 00:00:00 2001 From: rebeccamccabe Date: Sun, 1 Dec 2024 00:32:52 -0500 Subject: [PATCH 3/7] remove r_over_t since it is now unused --- mdocean/simulation/modules/geometry.m | 5 +---- mdocean/simulation/simulation.m | 4 ++-- 2 files changed, 3 insertions(+), 6 deletions(-) diff --git a/mdocean/simulation/modules/geometry.m b/mdocean/simulation/modules/geometry.m index 80098fca..af041e9e 100644 --- a/mdocean/simulation/modules/geometry.m +++ b/mdocean/simulation/modules/geometry.m @@ -1,5 +1,5 @@ function [V_d, m_m, m_f_tot, m_s_tot,... - A_c, A_lat_sub, r_over_t, ... + A_c, A_lat_sub, ... I, T, V_f_pct, V_s_pct, GM, mass,... CB_f_from_waterline,CG_f_from_waterline] = geometry(D_s, D_f, D_f_in, D_f_b, ... T_f_1, T_f_2, h_f, h_s, h_fs_clear, D_f_tu, ... @@ -144,9 +144,6 @@ A_c = [A_f_cross_top, A_vc_c, A_d_c]; A_lat_sub = [A_f_l A_vc_l A_d_l]; -r_over_t = [0,... % D_sft/(2*t_sf) - D_s/(2*t_s_r),... - 0];%D_or/(2*t_r)]; I = [I_f, I_vc, I_rp]; T = [T_f_2, T_s, h_d]; % drafts: used to calculated F_surge in dynamics.m m_m = m_f_m + m_s_m; % total mass of material diff --git a/mdocean/simulation/simulation.m b/mdocean/simulation/simulation.m index c5181d69..f3d6fe6d 100644 --- a/mdocean/simulation/simulation.m +++ b/mdocean/simulation/simulation.m @@ -4,7 +4,7 @@ %% Assemble inputs in = p; - + in.D_s = X(1); % inner diameter of float (m) in.D_f = X(2); % outer diameter of float (m) in.T_f_2 = X(3); % draft of float (m) @@ -32,7 +32,7 @@ %% Run modules [V_d, m_m, m_f_tot, m_s_tot, ... - A_c, A_lat_sub, r_over_t, ... + A_c, A_lat_sub, ... I, T, V_f_pct, V_s_pct, GM] = geometry(in.D_s, in.D_f, in.D_f_in, in.D_f_b, ... in.T_f_1, in.T_f_2, in.h_f, in.h_s, ... in.h_fs_clear, in.D_f_tu, in.t_f_t, ... From e74b22d149881f117e3bcc57ca0ce1be61c69c00 Mon Sep 17 00:00:00 2001 From: rebeccamccabe Date: Sun, 1 Dec 2024 04:58:59 -0500 Subject: [PATCH 4/7] add fatigue stress using operational sea states, separate from peak storm stress --- mdocean/inputs/parameters.m | 5 ++- mdocean/inputs/var_bounds.m | 6 ++- .../simulation/modules/dynamics/dynamics.m | 30 +++++++------ mdocean/simulation/modules/structures.m | 31 ++++++++++---- mdocean/simulation/simulation.m | 42 ++++++++++--------- 5 files changed, 69 insertions(+), 45 deletions(-) diff --git a/mdocean/inputs/parameters.m b/mdocean/inputs/parameters.m index 95eee7cb..a35e5b28 100644 --- a/mdocean/inputs/parameters.m +++ b/mdocean/inputs/parameters.m @@ -52,8 +52,9 @@ ... ...% Materials: [ Structural Steel ASTM-A36 (Ductile) ...% Concrete (Brittle) - ...% Stainless Steel 304 (Ductile) ] + ...% Stainless Steel 316 (Ductile) ASTM-A240 ] table("sigma_y","\sigma_y",{[36,4.5,30]* ksi2pa},"structures",true,"yield strength (Pa)"); + table("sigma_e","\sigma_e",{[58*.45, 0, 75*.45]*ksi2pa},"structures",true,"endurance limit (Pa)") table("rho_m","\rho_m",{[7850 2400 7900]},"structures",true,"material density (kg/m3)"); table("E","E",{[200e9, 5000*sqrt(4.5*ksi2pa), 200e9]},"structures",true,... "young's modulus (Pa)"); @@ -61,7 +62,7 @@ ...% RM3 CBS sheet 1.4 average of cells F21, F34, F46, F56 ...% https://www.concretenetwork.com/concrete-prices.html ...% https://agmetalminer.com/metal-prices/ - table("nu","\nu",0.3,"structures",true,"Poisson's ratio (-)"); + table("nu","\nu",{[0.36 0 0.29]},"structures",true,"Poisson's ratio (-)"); ... ...% Thicknesses and structures table("t_f_t","t_{f,t}",{0.50 * in2m},"structures",true,"float top thickness (m)"); diff --git a/mdocean/inputs/var_bounds.m b/mdocean/inputs/var_bounds.m index da6f7776..379e7fc9 100644 --- a/mdocean/inputs/var_bounds.m +++ b/mdocean/inputs/var_bounds.m @@ -83,7 +83,9 @@ b.X_start_struct = cell2struct(num2cell(b.X_starts),b.var_names(1:end-1)',1); b.constraint_names = {'float_too_heavy','float_too_light','spar_too_heavy','spar_too_light',... - 'stability','FOS_float_yield','FOS_col_yield','FOS_plate','FOS_col_buckling',... + 'stability','FOS_float_max','FOS_float_fatigue','FOS_col_max',... + 'FOS_col_fatigue','FOS_plate_max','FOS_plate_fatigue',... + 'FOS_col_local_max','FOS_col_local_fatigue',... 'pos_power','LCOE_max','irrelevant_max_force',... 'spar_height_up','spar_height_down','linear_theory'}; i1 = length(b.constraint_names); @@ -91,7 +93,7 @@ b.constraint_names{i} = strcat('prevent_slamming',num2str(i-i1)); end -b.lin_constraint_names = {'spar_natural_freq','float_spar_diam','float_spar_draft', +b.lin_constraint_names = {'spar_natural_freq','float_spar_diam','float_spar_draft',... 'float_spar_tops','float_seafloor','spar_seafloor'}; [~,idxs_sort] = sort(b.var_names(1:end-1)); % alphabetical design variable indices diff --git a/mdocean/simulation/modules/dynamics/dynamics.m b/mdocean/simulation/modules/dynamics/dynamics.m index 358c4e37..25432d8b 100644 --- a/mdocean/simulation/modules/dynamics/dynamics.m +++ b/mdocean/simulation/modules/dynamics/dynamics.m @@ -1,9 +1,12 @@ -function [F_heave_max, F_surge_max, F_ptrain_max, ... +function [F_heave_storm, F_surge_storm, F_heave_op, F_surge_op, F_ptrain_max, ... P_var, P_avg_elec, P_matrix_elec, X_constraints, B_p, X_u, X_f, X_s, P_matrix_mech] = dynamics(in,m_float,m_spar,V_d,draft) % use probabilistic sea states for power and PTO force and max amplitude [T,Hs] = meshgrid(in.T,in.Hs); - [P_matrix_mech,X_constraints,B_p,X_u,X_f,X_s,~,~,F_ptrain_max] = get_power_force(in,T,Hs,m_float,m_spar,V_d,draft); + [P_matrix_mech,X_constraints,B_p,... + X_u,X_f,X_s,F_heave_op,... + F_surge_op,F_ptrain_max] = get_power_force(in,T,Hs,m_float,m_spar,... + V_d,draft,in.F_max); % account for powertrain electrical losses P_matrix_elec = P_matrix_mech * in.eff_pto; @@ -18,8 +21,8 @@ assert(isreal(P_avg_elec)) % use max sea states for structural forces - [~,~,~,~,~,~,F_heave_max,F_surge_max,~] = get_power_force(in, ... - in.T_struct, in.Hs_struct, m_float, m_spar, V_d, draft); + [~,~,~,~,~,~,F_heave_storm,F_surge_storm,~] = get_power_force(in, ... + in.T_struct, in.Hs_struct, m_float, m_spar, V_d, draft, in.F_max); % coefficient of variance (normalized standard deviation) of power P_var = std(P_matrix_elec(:), in.JPD(:)) / P_avg_elec; @@ -29,7 +32,7 @@ end function [P_matrix, X_constraints, B_p, mag_X_u, mag_X_f, mag_X_s,... - F_heave_f, F_surge, F_ptrain_max] = get_power_force(in,T,Hs, m_float, m_spar, V_d, draft) + F_heave_f, F_surge, F_ptrain_max] = get_power_force(in,T,Hs, m_float, m_spar, V_d, draft, F_max) % get dynamic coefficients for float and spar % fixme: eventually should use in.D_f_in to allow a radial gap between float and spar @@ -52,7 +55,7 @@ mag_X_f,phase_X_f,... mag_X_s,phase_X_s,... B_p,K_p] = get_response_drag(w,m_f,m_s,m_c,B_h_f,B_h_s,B_c,K_h_f,K_h_s,... - F_f_mag,F_f_phase,F_s_mag,F_s_phase,in.F_max,... + F_f_mag,F_f_phase,F_s_mag,F_s_phase,F_max,... drag_const_f,drag_const_s,mag_v0_f,mag_v0_s, ... X_max,in.control_type,in.use_multibody,... in.X_tol,in.phase_X_tol,in.max_drag_iters); @@ -87,13 +90,14 @@ if nargout > 3 % powertrain force F_ptrain_max = max(mag_U,[],'all'); - F_ptrain_max = min(F_ptrain_max, in.F_max); + F_ptrain_max = min(F_ptrain_max, F_max); - % heave force: includes powertrain force and D'Alembert force - %F_heave_fund = sqrt( (mult .* B_p .* w).^2 + (mult .* K_p - m_float .* w.^2).^2 ) .* X; - - F_heave_f = combine_ptrain_dalembert_forces(m_float, w, mag_X_f, phase_X_f, mag_U, phase_U, in.F_max); - F_heave_s = combine_ptrain_dalembert_forces(m_spar, w, mag_X_s, phase_X_s, mag_U, phase_U, in.F_max); + % heave force: includes powertrain force and D'Alembert force + F_heave_f = combine_ptrain_dalembert_forces(m_float, w, mag_X_f, phase_X_f, mag_U, phase_U, F_max); + F_heave_s = combine_ptrain_dalembert_forces(m_spar, w, mag_X_s, phase_X_s, mag_U, phase_U, F_max); + + F_heave_f = max(F_heave_f,[],'all'); + F_heave_s = max(F_heave_s,[],'all'); % surge force F_surge = max(Hs,[],'all') * in.rho_w * in.g * V_d .* (1 - exp(-max(k_wvn,[],'all')*draft)); @@ -266,7 +270,7 @@ end % get force-saturated response - r = min(F_max ./ mag_U_unsat, 1);%fcn2optimexpr(@min, in.F_max ./ F_ptrain_unsat, 1); + r = min(F_max ./ mag_U_unsat, 1);%fcn2optimexpr(@min, F_max ./ F_ptrain_unsat, 1); alpha = 2/pi * ( 1./r .* asin(r) + sqrt(1 - r.^2) ); f_sat = alpha .* r; mult = get_multiplier(f_sat,m_f,B_f,K_f,w, B_f./B_p, K_f./K_p); % fixme this is wrong for multibody diff --git a/mdocean/simulation/modules/structures.m b/mdocean/simulation/modules/structures.m index e1d0d2df..83e8b8bc 100644 --- a/mdocean/simulation/modules/structures.m +++ b/mdocean/simulation/modules/structures.m @@ -1,6 +1,21 @@ function [FOS1Y,FOS2Y,FOS3Y,FOS_buckling] = structures(... - F_heave, F_surge, M, h_s, T_s, rho_w, g, ... - sigma_y, A_c, A_lat_sub, D_s, t_s_r, I, E, nu) + F_heave_storm, F_surge_storm, F_heave_op, F_surge_op, ... + M, h_s, T_s, rho_w, g, sigma_y, sigma_e, A_c, A_lat_sub, D_s, t_s_r, I, E, nu) + + F_heave_peak = max(F_heave_storm,F_heave_op); + F_surge_peak = max(F_surge_storm,F_surge_op); + + % peak + [FOS1Y, FOS2Y, FOS3Y, FOS_buckling] = structures_one_case(F_heave_peak, F_surge_peak, ... + M, h_s, T_s, rho_w, g, sigma_y(M), A_c, A_lat_sub, D_s, t_s_r, I, E, nu); + + % endurance limit (fatigue) + [FOS1Y(2), FOS2Y(2), FOS3Y(2), FOS_buckling(2)] = structures_one_case(F_heave_op, F_surge_op, ... + M, h_s, T_s, rho_w, g, sigma_e(M), A_c, A_lat_sub, D_s, t_s_r, I, E, nu); +end + +function [FOS1Y, FOS2Y, FOS3Y, FOS_spar_local] = structures_one_case(F_heave, F_surge, ... + M, h_s, T_s, rho_w, g, sigma_max, A_c, A_lat_sub, D_s, t_s_r, I, E, nu) %% Stress calculations depth = [0 T_s T_s]; % max depth @@ -27,17 +42,15 @@ sigma_vm = von_mises(sigma_rr, sigma_tt, sigma_zz, sigma_rt, sigma_tz, sigma_zr); %% Buckling calculation - FOS_spar = spar_combined_buckling(F_heave, E(M), I(2), h_s, D_s, A_c(2), t_s_r, ... - P_hydrostatic(2), sigma_y(M), nu); + [FOS_spar,FOS_spar_local] = spar_combined_buckling(F_heave, E(M), I(2), h_s, D_s, A_c(2), t_s_r, ... + P_hydrostatic(2), sigma_max, nu(M)); %% Factor of Safety (FOS) Calculations - FOS_yield = sigma_y(M) ./ sigma_vm; + FOS_yield = sigma_max ./ sigma_vm; FOS1Y = FOS_yield(1); FOS2Y = FOS_spar; FOS3Y = FOS_yield(3); - FOS_buckling = 3; - end function [FOS_spar, FOS_spar_local] = spar_combined_buckling(F, E, I, L, D, A, t, q, sigma_0, nu) @@ -51,7 +64,7 @@ % local buckling of plate element: 2/9.7 k_s = 1.33; % uniform compression, fixed-free, from table 3 on page 30 - s = D; % not sure if this is correct + s = D; % fixme - not sure if this is correct P_r = 0.6; % steel proportional linear elastic limit sigma_Ex = k_s * pi^2 * E / (12 * (1-nu^2)) * (t/s)^2; if sigma_Ex <= P_r * sigma_0 @@ -85,7 +98,7 @@ % local buckling of tube: final part of 2/7.3, which combines 2/9.1 % and 2/9.5 - FOS_spar_local = 0; + FOS_spar_local = 3; % fixme: need to implement end function s_vm = von_mises(s_11, s_22, s_33, s_12, s_23, s_31) diff --git a/mdocean/simulation/simulation.m b/mdocean/simulation/simulation.m index f3d6fe6d..b28d04e9 100644 --- a/mdocean/simulation/simulation.m +++ b/mdocean/simulation/simulation.m @@ -43,40 +43,44 @@ m_f_tot = max(m_f_tot,1e-3); % zero out negative mass produced by infeasible inputs -[F_heave_max, F_surge_max, F_ptrain_max, ... - P_var, P_avg_elec, P_matrix_elec, ... - X_constraints] = dynamics(in, m_f_tot, m_s_tot, V_d, T); +[F_heave_storm, F_surge_storm, ... + F_heave_op, F_surge_op, F_ptrain_max, ... + P_var, P_avg_elec, P_matrix_elec, ... + X_constraints] = dynamics(in, m_f_tot, m_s_tot, V_d, T); [FOS_float,FOS_spar,FOS_damping_plate,... - FOS_spar_local_buckling] = structures(... - F_heave_max, F_surge_max,... - in.M, in.h_s, in.T_s, in.rho_w, in.g, in.sigma_y, A_c, ... - A_lat_sub, in.D_s, in.t_s_r, I, in.E, in.nu); + FOS_spar_local] = structures(F_heave_storm, F_surge_storm, F_heave_op, F_surge_op, ... + in.M, in.h_s, in.T_s, in.rho_w, in.g, in.sigma_y, in.sigma_e, ... + A_c, A_lat_sub, in.D_s, in.t_s_r, I, in.E, in.nu); LCOE = econ(m_m, in.M, in.cost_m, in.N_WEC, P_avg_elec, in.FCR, in.cost_perN, in.F_max, in.eff_array); %% Assemble constraints g(x) >= 0 -num_g = 15+numel(p.JPD); +num_g = 19+numel(p.JPD); g = zeros(1,num_g); g(1) = V_f_pct; % prevent float too heavy g(2) = 1 - V_f_pct; % prevent float too light g(3) = V_s_pct; % prevent spar too heavy g(4) = 1 - V_s_pct; % prevent spar too light g(5) = GM; % pitch stability of float-spar system -g(6) = FOS_float / p.FOS_min - 1; % float survives max force -g(7) = FOS_spar / p.FOS_min - 1; % spar survives max force -g(8) = FOS_damping_plate / p.FOS_min - 1; % damping plate survives max force -g(9) = FOS_spar_local_buckling / p.FOS_min - 1; % spar survives max force in buckling -g(10) = P_avg_elec; % positive power +g(6) = FOS_float(1) / p.FOS_min - 1; % float survives max force +g(7) = FOS_float(2) / p.FOS_min - 1; % float survives fatigue +g(8) = FOS_spar(1) / p.FOS_min - 1; % spar survives max force +g(9) = FOS_spar(2) / p.FOS_min - 1; % spar survives fatigue +g(10) = FOS_damping_plate(1) / p.FOS_min - 1; % damping plate survives max force +g(11) = FOS_damping_plate(2) / p.FOS_min - 1; % damping plate survives fatigue +g(12) = FOS_spar_local(1) / p.FOS_min - 1; % spar survives max force in local buckling +g(13) = FOS_spar_local(2) / p.FOS_min - 1; % spar survives fatigue in local buckling +g(14) = P_avg_elec; % positive power %1 + min(Kp_over_Ks,[],'all'); % spar heave stability (positive effective stiffness) -g(11) = p.LCOE_max/LCOE - 1; % prevent more expensive than threshold -g(12) = F_ptrain_max/in.F_max - 1; % prevent irrelevant max force - +g(15) = p.LCOE_max/LCOE - 1; % prevent more expensive than threshold +g(16) = F_ptrain_max/in.F_max - 1; % prevent irrelevant max force - % this constraint should always be active % and is only required when p.cost_perN = 0. -g(13) = X_constraints(1); % prevent float rising above top of spar -g(14) = X_constraints(2); % prevent float going below bottom of spar -g(15) = X_constraints(3); % float amplitude obeys linear theory -g(16:end) = X_constraints(4:end); % prevent rising out of water/slamming +g(17) = X_constraints(1); % prevent float rising above top of spar +g(18) = X_constraints(2); % prevent float going below bottom of spar +g(19) = X_constraints(3); % float amplitude obeys linear theory +g(20:end) = X_constraints(4:end); % prevent rising out of water/slamming % fixme: add another X_constraint that h_fs_clear is greater than X_u criteria = all(~isinf(g)) && all(~isnan(g)) && all(isreal(g)); From 7a7f33b7e89ca7ed91d0bf5cdb089bea35ac33d5 Mon Sep 17 00:00:00 2001 From: rebeccamccabe Date: Sun, 1 Dec 2024 05:06:03 -0500 Subject: [PATCH 5/7] fix renamed variable --- mdocean/simulation/simulation.m | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mdocean/simulation/simulation.m b/mdocean/simulation/simulation.m index b28d04e9..f59e477f 100644 --- a/mdocean/simulation/simulation.m +++ b/mdocean/simulation/simulation.m @@ -109,7 +109,7 @@ val.LCOE = LCOE; val.power_avg = P_avg_elec; val.power_max = max(P_matrix_elec,[],'all'); - val.force_heave = F_heave_max; + val.force_heave = F_heave_storm; val.force_ptrain = F_ptrain_max; val.FOS_spar = FOS_spar; val.c_v = P_var; From 67a06d9692c25eaad5a31a3249ac7c21259574b3 Mon Sep 17 00:00:00 2001 From: rebeccamccabe Date: Sun, 1 Dec 2024 05:16:28 -0500 Subject: [PATCH 6/7] make val.FOS_spar scalar again for validation to work --- mdocean/simulation/simulation.m | 2 +- tests/test.m | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/mdocean/simulation/simulation.m b/mdocean/simulation/simulation.m index f59e477f..079fe4df 100644 --- a/mdocean/simulation/simulation.m +++ b/mdocean/simulation/simulation.m @@ -111,7 +111,7 @@ val.power_max = max(P_matrix_elec,[],'all'); val.force_heave = F_heave_storm; val.force_ptrain = F_ptrain_max; - val.FOS_spar = FOS_spar; + val.FOS_spar = FOS_spar(1); val.c_v = P_var; val.B_p = B_p; val.X_u = X_u; diff --git a/tests/test.m b/tests/test.m index 3b3163ee..7e73c13c 100644 --- a/tests/test.m +++ b/tests/test.m @@ -4,7 +4,7 @@ % class based unit tests, as in https://www.mathworks.com/help/matlab/matlab_prog/class-based-unit-tests.html properties (Constant) - run_slow_tests = true; + run_slow_tests = false; slow_figs = [6 7 8]; slow_tabs = 7; end From 1e68b64e06e31442b56e9b8cfe400ce21b5d3d79 Mon Sep 17 00:00:00 2001 From: rebeccamccabe Date: Sun, 1 Dec 2024 05:21:00 -0500 Subject: [PATCH 7/7] turn slow tests back on (got turned off accidentally) --- tests/test.m | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test.m b/tests/test.m index 7e73c13c..3b3163ee 100644 --- a/tests/test.m +++ b/tests/test.m @@ -4,7 +4,7 @@ % class based unit tests, as in https://www.mathworks.com/help/matlab/matlab_prog/class-based-unit-tests.html properties (Constant) - run_slow_tests = false; + run_slow_tests = true; slow_figs = [6 7 8]; slow_tabs = 7; end