Merge pull request #106 from byuflowlab/optimization

minor bug fixes
byuflowlab · Mar 29, 2023 · 32456d5 · 32456d5
2 parents 6881dbf + f441746
commit 32456d5
Show file tree

Hide file tree

Showing 4 changed files with 8 additions and 23 deletions.
diff --git a/docs/src/literate/dynamic-joined-wing.jl b/docs/src/literate/dynamic-joined-wing.jl
@@ -217,19 +217,6 @@ end
 #md # ![](../assets/dynamic-joined-wing-M2.svg)
 #md # ![](../assets/dynamic-joined-wing-M3.svg)
 
-#md # ![](../assets/dynamic-joined-wing-u1.svg)
-#md # ![](../assets/dynamic-joined-wing-u2.svg)
-#md # ![](../assets/dynamic-joined-wing-u3.svg)
-#md # ![](../assets/dynamic-joined-wing-theta1.svg)
-#md # ![](../assets/dynamic-joined-wing-theta2.svg)
-#md # ![](../assets/dynamic-joined-wing-theta3.svg)
-#md # ![](../assets/dynamic-joined-wing-F1.svg)
-#md # ![](../assets/dynamic-joined-wing-F2.svg)
-#md # ![](../assets/dynamic-joined-wing-F3.svg)
-#md # ![](../assets/dynamic-joined-wing-M1.svg)
-#md # ![](../assets/dynamic-joined-wing-M2.svg)
-#md # ![](../assets/dynamic-joined-wing-M3.svg)
-
 #-
 
 # These graphs are identical to those presented in "GEBT: A general-purpose nonlinear analysis tool for composite beams" by Wenbin Yu and Maxwell Blair.

diff --git a/docs/src/literate/rotating.jl b/docs/src/literate/rotating.jl
@@ -129,6 +129,7 @@ plot(
     xlim = (0, 760),
     xticks = 0:100:750,
     xlabel = "Angular Speed (RPM)",
+    ylim = (0, 12)
     yticks = 0.0:2:12,
     ylabel = "\$M_z\$ at the root (lb-in)",
     grid = false,

diff --git a/src/analyses.jl b/src/analyses.jl
@@ -642,7 +642,7 @@ function steady_output!(system, x, p, constants)
     assembly, pcond, dload, pmass, gvec, vb_p, ωb_p, ab_p, αb_p = steady_parameters(x, p, constants)
 
     # initialize state rate vector (if necessary)
-    dx = typeof(system.dx) <: typeof(x) ? system.dx .= 0 : similar(x) .= system.x
+    dx = typeof(system.dx) <: typeof(x) ? system.dx .= 0 : zero(x)
 
     # extract body frame accelerations
     ab, αb = body_accelerations(x, indices.icol_body, ab_p, αb_p)
@@ -1002,7 +1002,7 @@ function left_eigenvectors(K, M, λ, V)
     for iλ = 1:nev
 
         # factorize (K + λ*M)'
-        KmλMfact = factorize(K' + λ[iλ]'*M')
+        KmλMfact = factorize(K' + (λ[iλ]+1e-9)'*M')
 
         # initialize left eigenvector
         for i = 1:nx
@@ -1050,7 +1050,7 @@ function left_eigenvectors(K, M::SparseMatrixCSC, λ, V)
     for iλ = 1:nev
 
         # factorize (K + λ*M)'
-        KmλMfact = factorize(K' + λ[iλ]'*M')
+        KmλMfact = factorize(K' + (λ[iλ]+1e-9)'*M')
 
         # initialize left eigenvector
         for i = 1:nx
@@ -2353,7 +2353,7 @@ function time_domain_analysis!(system::DynamicSystem, assembly, tvec;
     # --- Initialize Time-Domain Solution --- #
 
     # initialize storage for each time step
-    history = Vector{AssemblyState{eltype(x)}}(undef, length(save))
+    history = Vector{AssemblyState{eltype(x), Vector{PointState{eltype(x)}}, Vector{ElementState{eltype(x)}}}}(undef, length(save))
     isave = 1
 
     # add initial state to the solution history
@@ -2496,7 +2496,7 @@ function newmark_parameters(p, constants)
     Ωdot_init = [SVector{3}(p_i[12*(ip-1)+10], p_i[12*(ip-1)+11], p_i[12*(ip-1)+12]) for ip = 1:length(assembly.points)]
 
     # overwrite default assembly and parameters (if applicable)
-    parameters = isnothing(xpfunc) ? pfunc(p, t) : xpfunc(x, p, t)
+    parameters = isnothing(xpfunc) ? pfunc(p_p, t) : xpfunc(x, p_p, t)
     assembly = get(parameters, :assembly, assembly)
     prescribed_conditions = get(parameters, :prescribed_conditions, prescribed_conditions)
     distributed_loads = get(parameters, :distributed_loads, distributed_loads)
@@ -2605,9 +2605,6 @@ function newmark_output(system, x, p, constants)
     # combine constants and parameters
     assembly, pcond, dload, pmass, gvec, vb_p, ωb_p, udot_init, θdot_init, Vdot_init, Ωdot_init, dt = newmark_parameters(p, constants)
 
-    # update acceleration state variable indices
-    update_body_acceleration_indices!(indices, pcond)
-
     # initialize state rate vector (if necessary)
     dx = typeof(system.dx) <: typeof(x) ? system.dx .= 0 : similar(x) .= 0
 

diff --git a/src/output.jl b/src/output.jl
@@ -89,7 +89,7 @@ Base.eltype(::Type{AssemblyState{TF, TP, TE}}) where {TF, TP, TE} = TF
 function AssemblyState{TF, TP, TE}(a::AssemblyState) where {TF, TP, TE}
     return AssemblyState{TF, TP, TE}(a.points, a.elements)
 end
-Base.convert(::Type{AssemblyState{TF, TP, TE}}, p::AssemblyState) where {TF, TP, TE} = AssemblyState{TF, TP, TE}(p)
+Base.convert(::Type{AssemblyState{TF, TP, TE}}, a::AssemblyState) where {TF, TP, TE} = AssemblyState{TF, TP, TE}(a)
 
 """
     AssemblyState(system, assembly; prescribed_conditions = Dict())
@@ -445,7 +445,7 @@ function extract_element_state(dx, x, system::ExpandedSystem, assembly, ielem;
     V = CtCab*V
     Ω = CtCab*Ω
 
-    # rotate linear and angular veocity rates into the body frame
+    # rotate linear and angular velocity rates into the body frame
     CtCabdot = tilde(C'*get_Q(θ)*θdot)*CtCab
     Vdot = CtCabdot*V + CtCab*Vdot
     Ωdot = CtCabdot*Ω + CtCab*Ωdot