Merge pull request #42 from sandialabs/eparish1/opinf-add-sidesets

added code for multiple BCs per sideset

examples/ahead/overlap/cuboid/dynamic-opinf-fom/make_op_inf_models.py

@@ -20,13 +20,13 @@
     displacement_snapshots[free_dofs[:,:]==False] = 0.
 
     #Get sideset snapshots
-    sidesets = ["nsx-","nsy-","ssz-","nsz+"]
+    sidesets = ["nsx--x","nsy--y","ssz-","nsz+-z"]
     sideset_snapshots = normaopinf.readers.load_sideset_displacement_csv_files(solution_directory=cur_dir,sidesets=sidesets,solution_id=solution_id,skip_files=1)
 
-
+    tolerance = 1.e-5
     # Create bases for sidesets 
-    my_energy_truncater = romtools.vector_space.utils.EnergyBasedTruncater(1. - 1.e-5)
-    
+    my_energy_truncater = romtools.vector_space.utils.EnergyBasedTruncater(1. - tolerance)
+
     ss_tspace = {}
     reduced_sideset_snapshots = {}
     for sideset in sidesets:
@@ -80,7 +80,7 @@
     uhat = romtools.rom.optimal_l2_projection(displacement_snapshots,trial_space)
     u_ddots = normaopinf.calculus.d2dx2(displacement_snapshots,times)
     uhat_ddots = romtools.rom.optimal_l2_projection(u_ddots*1.,trial_space)
-    l2solver = opinf.lstsq.L2Solver(regularizer=1e-8)
+    l2solver = opinf.lstsq.L2Solver(regularizer=5e-3)
     opinf_model = opinf.models.ContinuousModel("AB",solver=l2solver)
 
     opinf_model.fit(states=uhat, ddts=uhat_ddots,inputs=reduced_stacked_sideset_snapshots)

src/ics_bcs.jl

@@ -221,7 +221,19 @@ function apply_bc(model::LinearOpInfRom, bc::SMDirichletBC)
         disp_val = model.fom_model.current[bc.offset,node_index] - model.fom_model.reference[bc.offset, node_index]
         push!(bc_vector,disp_val)
     end
-    op_name = "B_"*bc.node_set_name 
+
+    offset = bc.offset
+    if offset == 1
+      offset_name = "x"
+    end
+    if offset == 2
+      offset_name = "y"
+    end
+    if offset == 3
+      offset_name = "z"
+    end
+
+    op_name = "B_"*bc.node_set_name  * "-" * offset_name
     bc_operator = model.opinf_rom[op_name] 
     # SM Dirichlet BC are only defined on a single x,y,z
     model.reduced_boundary_forcing[:] += bc_operator[1,:,:] * bc_vector

src/time_integrator.jl

@@ -509,6 +509,8 @@ function write_step_csv(integrator::TimeIntegrator, model::SolidMechanics, sim_i
         refe_filename = sim_id_string * "refe" * ".csv"
         writedlm_nodal_array(refe_filename, model.reference)
     end
+    free_dofs_filename = sim_id_string * "free_dofs" * index_string * ".csv"
+    writedlm(free_dofs_filename, model.free_dofs)
     curr_filename = sim_id_string * "curr" * index_string * ".csv"
     writedlm_nodal_array(curr_filename, model.current)
     disp_filename = sim_id_string * "disp" * index_string * ".csv"
@@ -527,7 +529,7 @@ function write_step_csv(integrator::TimeIntegrator, model::SolidMechanics, sim_i
     end
 end
 
-function write_sideset_step_csv(params::Dict{Any,Any},integrator::DynamicTimeIntegrator, model::SolidMechanics, sim_id::Integer)
+function write_sideset_step_csv(params::Dict{String,Any},integrator::DynamicTimeIntegrator, model::SolidMechanics, sim_id::Integer)
   stop = integrator.stop
   index_string = "-" * string(stop, pad = 4)
   sim_id_string = string(sim_id, pad = 2) * "-"
@@ -537,20 +539,30 @@ function write_sideset_step_csv(params::Dict{Any,Any},integrator::DynamicTimeInt
   for bc ∈ model.boundary_conditions
       if (typeof(bc) == SMDirichletBC)
         node_set_name = bc.node_set_name 
-        curr_filename = sim_id_string * node_set_name * "curr" * index_string * ".csv"
-        disp_filename = sim_id_string * node_set_name * "disp" * index_string * ".csv"
-        velo_filename = sim_id_string * node_set_name * "velo" * index_string * ".csv"
-        acce_filename = sim_id_string * node_set_name * "acce" * index_string * ".csv"
+        offset = bc.offset
+        if offset == 1
+          offset_name = "x"
+        end
+        if offset == 2
+          offset_name = "y"
+        end
+        if offset == 3
+          offset_name = "z"
+        end
+        curr_filename = sim_id_string * node_set_name * "-" * offset_name * "-curr" * index_string * ".csv"
+        disp_filename = sim_id_string * node_set_name * "-" * offset_name * "-disp" * index_string * ".csv"
+        velo_filename = sim_id_string * node_set_name * "-" * offset_name * "-velo" * index_string * ".csv"
+        acce_filename = sim_id_string * node_set_name * "-" * offset_name * "-acce" * index_string * ".csv"
         writedlm(curr_filename, model.current[bc.offset,bc.node_set_node_indices])
         writedlm(velo_filename, model.velocity[bc.offset,bc.node_set_node_indices])
         writedlm(acce_filename, model.acceleration[bc.offset,bc.node_set_node_indices])
         writedlm(disp_filename, model.current[bc.offset,bc.node_set_node_indices] - model.reference[bc.offset,bc.node_set_node_indices])
       elseif (typeof(bc) == SMOverlapSchwarzBC)
         side_set_name = bc.side_set_name 
-        curr_filename = sim_id_string * side_set_name * "curr" * index_string * ".csv"
-        disp_filename = sim_id_string * side_set_name * "disp" * index_string * ".csv"
-        velo_filename = sim_id_string * side_set_name * "velo" * index_string * ".csv"
-        acce_filename = sim_id_string * side_set_name * "acce" * index_string * ".csv"
+        curr_filename = sim_id_string * side_set_name * "-curr" * index_string * ".csv"
+        disp_filename = sim_id_string * side_set_name * "-disp" * index_string * ".csv"
+        velo_filename = sim_id_string * side_set_name * "-velo" * index_string * ".csv"
+        acce_filename = sim_id_string * side_set_name * "-acce" * index_string * ".csv"
         writedlm_nodal_array(curr_filename, model.current[:,bc.side_set_node_indices])
         writedlm_nodal_array(velo_filename, model.velocity[:,bc.side_set_node_indices])
         writedlm_nodal_array(acce_filename, model.acceleration[:,bc.side_set_node_indices])
@@ -559,6 +571,7 @@ function write_sideset_step_csv(params::Dict{Any,Any},integrator::DynamicTimeInt
   end
 end
 
+
 function write_step_csv(integrator::DynamicTimeIntegrator, model::OpInfModel, sim_id::Integer)
     stop = integrator.stop
     index_string = "-" * string(stop, pad = 4)

test/opinf-schwarz-overlap-cuboid-hex8.jl

@@ -23,7 +23,7 @@
     rm("cuboid-2.e")
     rm("opinf-operator.npz")
     print(model_coarse.reduced_state[1])
-    @test model_coarse.reduced_state[1] ≈ 0.04374576440784916 rtol = 1.0e-06
-    @test model_coarse.reduced_state[2] ≈ 0.04374576440786179 rtol = 1.0e-06
-    @test model_coarse.reduced_state[3] ≈ 0.19557108653526997 rtol = 1.0e-06
+    @test model_coarse.reduced_state[1] ≈ 0.04508462384388636 rtol = 1.0e-06
+    @test model_coarse.reduced_state[2] ≈ 0.0450846238438998 rtol = 1.0e-06
+    @test model_coarse.reduced_state[3] ≈ 0.21766337487648266 rtol = 1.0e-06
 end