Skip to content

Commit

Permalink
Merge cda03bc into 87558c7
Browse files Browse the repository at this point in the history
  • Loading branch information
relf committed Apr 8, 2022
2 parents 87558c7 + cda03bc commit 2d2d251
Show file tree
Hide file tree
Showing 12 changed files with 217 additions and 42 deletions.
4 changes: 2 additions & 2 deletions .github/workflows/tests.yml
Original file line number Diff line number Diff line change
Expand Up @@ -42,10 +42,10 @@ jobs:
- name: Make Thrift Compiler
if: steps.cache-thrift.outputs.cache-hit != 'true'
run: |
export THRIFT_VERSION=0.15.0
export THRIFT_VERSION=0.16.0
export buildDeps="automake bison curl flex g++ libboost-dev libboost-filesystem-dev libboost-program-options-dev libboost-system-dev libboost-test-dev libevent-dev libssl-dev libtool make pkg-config"
sudo apt-get install -y --no-install-recommends $buildDeps && sudo rm -rf /var/lib/apt/lists/*
curl -sSL "https://dlcdn.apache.org/thrift/$THRIFT_VERSION/thrift-$THRIFT_VERSION.tar.gz" -o thrift.tar.gz
curl -sSL "http://archive.apache.org/dist/thrift/$THRIFT_VERSION/thrift-$THRIFT_VERSION.tar.gz" -o thrift.tar.gz
mkdir -p thrift
mkdir -p install
tar zxf thrift.tar.gz -C ${PWD}/thrift --strip-components=1
Expand Down
3 changes: 3 additions & 0 deletions .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -73,3 +73,6 @@ services/whatsopt_server/optimizer_store/oneramdao

# Ignore XDSMjs
app/javascript/XDSMjs

# Ignore Openmdao reports
**/reports
1 change: 1 addition & 0 deletions app/lib/whats_opt/egmdo_generator.rb
Original file line number Diff line number Diff line change
Expand Up @@ -28,6 +28,7 @@ def _generate_code(gendir, options = {})
_generate("doe_factory.py", "egmdo/doe_factory.py.erb", egmdo_dir)
_generate("gp_factory.py", "egmdo/gp_factory.py.erb", egmdo_dir)
_generate("random_analysis.py", "egmdo/random_analysis.py.erb", egmdo_dir)
_generate("random_vec_analysis.py", "egmdo/random_vec_analysis.py.erb", egmdo_dir)
_generate("__init__.py", "__init__.py.erb", egmdo_dir)
if @driver_name
@driver = OpenmdaoDriverFactory.new(@driver_name, @driver_options).create_driver
Expand Down
35 changes: 16 additions & 19 deletions app/lib/whats_opt/templates/egmdo/algorithms.py.erb
Original file line number Diff line number Diff line change
Expand Up @@ -10,26 +10,21 @@ from openmdao_extensions.smt_doe_driver import SmtDOEDriver
from openmdao_extensions.openturns_doe_driver import OpenturnsDOEDriver
from openmdao.utils.mpi import MPI

from <%= @pkg_prefix %>egmdo.random_analysis import <%= @mda.py_classname %>RandomAnalysis
from <%= @pkg_prefix %>egmdo.random_vec_analysis import <%= @mda.py_classname %>RandomVecAnalysis

from <%= @pkg_prefix %>egmdo.doe_factory import DoeFactory
from <%= @pkg_prefix %>egmdo.gp_factory import GpFactory

def compute_doe(discipline_factory, gp_factory, design_vars, n_cases, out_sqlitefilename, parallel=False):
pb = om.Problem(<%= @mda.py_classname %>RandomAnalysis(discipline_factory, gp_factory))
dists = []
<% dim = 0 %>
<%- @mda.egmdo_random_variables.each do |param| -%>
# _xi_<%= param.py_varname %>
dists.append(ot.Normal(0.0, 1.0))
<%- dim = dim + 1 -%>
<%- end -%>
pb = om.Problem(<%= @mda.py_classname %>RandomVecAnalysis(discipline_factory, gp_factory, n_cases))

dists = [ot.Normal(0.0,1.0)] * n_cases * <%= @mda.egmdo_random_variables.size %>

# Dependency between variables can be specified by choosing a specific copula
copula = ot.IndependentCopula(<%= dim %>) # default to no dependency
copula = ot.IndependentCopula(n_cases * <%= @mda.egmdo_random_variables.size %>) # default to no dependency

comp_dist = ot.ComposedDistribution(dists, copula)
pb.driver = OpenturnsDOEDriver(n_samples=n_cases, distribution=comp_dist)
pb.driver = OpenturnsDOEDriver(n_samples=1, distribution=comp_dist)
pb.driver.options['run_parallel'] = parallel

if parallel and MPI:
Expand All @@ -55,7 +50,7 @@ def compute_doe(discipline_factory, gp_factory, design_vars, n_cases, out_sqlite
for name, val in design_vars.items():
pb[name] = val
<%- @mda.egmdo_random_variables.each do |param| -%>
pb['_xi_<%= param.py_varname %>'] = 0.0
pb['_xi_<%= param.py_varname %>'] = np.zeros(n_cases)
<%- end -%>
pb.run_driver()
pb.cleanup()
Expand All @@ -76,27 +71,29 @@ def read_doe(qoi_names, n_doe_pce, outfile, parallel):
doe = np.zeros((n_doe_pce, n_des_var + len(qoi_names)))
coupling_var_doe = np.zeros((n_doe_pce, n_des_var))

obj_name = qoi_names[0]
obj_name = qoi_names[-1]
offset = 0
for file in files:
cr = om.CaseReader(file)
driver_cases = cr.list_cases("driver", out_stream=None)
case = cr.get_case(driver_cases[0])
n_doe = len(driver_cases)
print(f"{file} {n_doe}")
if n_doe == 0:
continue
case = cr.get_case(driver_cases[0])

for k in range(n_doe):
i = offset + k
case = cr.get_case(driver_cases[k])
for j in range(n_des_var):
doe[i, j] = case.outputs[des_var_names[j]]
doe[i, -1] = case[obj_name]
doe[:n_doe_pce, j] = case.outputs[des_var_names[j]]
doe[:n_doe_pce, -1] = case[obj_name][:, 0]
<% nb_cstrs = @mda.constraint_variables.size %>
<%- @mda.constraint_variables.each_with_index do |cstr_var, k| -%>
doe[i, <%= -(nb_cstrs+1)+k %>] = case["<%= cstr_var.py_varname %>"]
doe[:n_doe_pce, <%= -(nb_cstrs+1)+k %>] = case["<%= cstr_var.py_varname %>"][:,0]
<%- end -%>
<%- @mda.coupling_variables.each_with_index do |cv, j| -%>
coupling_var_doe[i, <%= j %>] = case["<%= cv.py_varname %>"]
coupling_var_doe[:n_doe_pce, <%= j %>] = case["<%= cv.py_varname %>"][:,0]
<%- end -%>

offset += n_doe
Expand All @@ -105,7 +102,7 @@ def read_doe(qoi_names, n_doe_pce, outfile, parallel):


def analyze_sensitivity(qoi_names, doe, threshold_coeff_var=1e-3):
obj_name = qoi_names[0]
obj_name = qoi_names[-1]

print(f"Estimated mean of '{obj_name}' = {doe[:, -1].mean()}", )
print(f"Estimated coeff of variation of '{obj_name}' = {doe[:, -1].std()/doe[:, -1].mean()}")
Expand Down
7 changes: 4 additions & 3 deletions app/lib/whats_opt/templates/egmdo/doe_factory.py.erb
Original file line number Diff line number Diff line change
Expand Up @@ -121,7 +121,8 @@ class DoeFactory:
for i in range(n_design_vars):
dim_design_vars.append(len(case.get_design_vars()[design_vars_names[i]]))

problem = om.Problem(self.discipline_factory.create_<%= disc.snake_modulename %>(), comm=self.comm)
problem = om.Problem(comm=self.comm)
problem.model.add_subsystem("<%= disc.py_classname %>", self.discipline_factory.create_<%= disc.snake_modulename %>(), promotes=["*"])
problem.setup()
<%- disc.design_variables.each do |var| -%>
problem['<%= var.py_varname %>'] = design_vars['<%= var.py_varname %>']
Expand All @@ -133,12 +134,12 @@ class DoeFactory:
# run the analysis
problem.run_model()

list_inputs = problem.model.list_inputs(out_stream=None)
list_inputs = problem.model.list_inputs(out_stream=None, prom_name=True)
# update the DoE with the new inputs and update the inputs of the problem with
# the correct value of the coupling variables

for inp in list_inputs:
name = inp[0]
name = inp[1]['prom_name']
ind = design_vars_names.index(name)
start = int(np.array(dim_design_vars[:ind + 1]).sum())
# we look for inp in the dict coupling_vars
Expand Down
13 changes: 6 additions & 7 deletions app/lib/whats_opt/templates/egmdo/gp_factory.py.erb
Original file line number Diff line number Diff line change
Expand Up @@ -20,13 +20,13 @@ class GpFactory:
# load or create the DoE
doe_file = self.doe_factory.doe_numpy_filename("<%= disc.snake_modulename %>")
if (not os.path.exists(doe_file)):
print(f"Create DOE for <%= disc.snake_modulename %> with {n_cases} samples: {doe_file}")
print(f"Create DOE for '<%= disc.name %>' with {n_cases} samples: {doe_file}")
self.doe_factory.create_doe_<%= disc.snake_modulename %>(n_cases)
doe = np.load(doe_file)

# create the gp surrogate model
gp_file = self.gp_filename("<%= disc.snake_modulename %>", "<%= v.py_varname %>")
print(f"Build GP for <%= disc.snake_modulename %> : {gp_file}")
print(f"Build GP for '<%= disc.name %>' : {gp_file}")
input_dim=<%= disc.input_variables.map(&:dim).inject(0, :+) %>
gp = KRG(theta0=[1e-2] * input_dim, print_global=False)
gp.set_training_values(doe[:, :-<%= outdim %>], doe[:, -<%= outdim-i %>])
Expand All @@ -35,17 +35,16 @@ class GpFactory:
pickle.dump(gp, f)

def update_gp_<%= disc.snake_modulename %>_<%= v.py_varname %>(self, design_vars, coupling_vars):
print(f"Update DOE for <%= disc.snake_modulename %>")
print(f"Update DOE for '<%= disc.name %>'")
new_doe = self.doe_factory.update_doe_<%= disc.snake_modulename %>(design_vars, coupling_vars)
self.create_gp_<%= disc.snake_modulename %>_<%= v.py_varname %>()
return new_doe
<%- end -%>
<%- end -%>

def create_all_gps(self, n_cases=4):
<%- @mda.egmdo_random_variables.each do |v| -%>
self.create_gp_<%= v.discipline.snake_modulename %>_<%= v.py_varname %>(n_cases)
<%- @mda.egmdo_random_disciplines.each do |disc| -%>
<%- disc.output_variables.each do |v| -%>
self.create_gp_<%= disc.snake_modulename %>_<%= v.py_varname %>(n_cases)
<%- end -%>
<%- if @mda.egmdo_random_disciplines.empty? -%>
pass
<%- end -%>
8 changes: 5 additions & 3 deletions app/lib/whats_opt/templates/egmdo/openmdao_egmda.py.erb
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ from numpy import nan, inf
import openmdao.api as om
from mda_init import initialize
from <%= @pkg_prefix %><%= @mda.py_modulename %> import <%= @mda.py_classname %>Factory as DisciplineFactory
from <%= @pkg_prefix %>egmdo.random_vec_analysis import <%= @mda.py_classname %>RandomVecAnalysis
from <%= @pkg_prefix %>egmdo.random_analysis import <%= @mda.py_classname %>RandomAnalysis
from <%= @pkg_prefix %>egmdo.doe_factory import DoeFactory
from <%= @pkg_prefix %>egmdo.gp_factory import GpFactory
Expand Down Expand Up @@ -57,8 +58,9 @@ class <%= @mda.py_classname %>EgDiscipline(om.ExplicitComponent):
parallel=self.parallel
)

# After the convergence we run the problem at the mean value of the disciplinary surrogate models to fill the outputs
pb = om.Problem(<%= @mda.py_classname %>RandomAnalysis(self.discipline_factory, self.gp_factory))
# After the convergence we run the problem at the mean value of the disciplinary
# surrogate models to fill the outputs
pb = om.Problem(<%= @mda.py_classname %>RandomVecAnalysis(self.discipline_factory, self.gp_factory))
# pb.model.nonlinear_solver.options['err_on_non_converge'] = False
# pb.model.linear_solver.options['err_on_non_converge'] = False
pb.model.nonlinear_solver.options['iprint'] = 0
Expand All @@ -69,7 +71,7 @@ class <%= @mda.py_classname %>EgDiscipline(om.ExplicitComponent):
<% end %>
pb.run_model()
<% @mda.output_variables.each do |var| -%>
outputs['<%= var.name %>'] = pb['<%= var.name %>']
outputs['<%= var.name %>'] = pb['<%= var.name %>'].mean()
<% end %>

def compute_partials(self, inputs, partials):
Expand Down
4 changes: 2 additions & 2 deletions app/lib/whats_opt/templates/egmdo/random_analysis.py.erb
Original file line number Diff line number Diff line change
Expand Up @@ -47,7 +47,7 @@ class <%= disc.py_classname %>RandomDiscipline(om.ExplicitComponent):
<% end %>
for name in sorted(inputs.keys()):
if name not in random_inputs:
inputs_gp = np.concatenate((inputs_gp, inputs[name]))
inputs_gp = np.concatenate((inputs_gp, inputs[name].flatten()))
inputs_gp = np.atleast_2d(inputs_gp)
<% disc.output_variables.each do |v|%>
sigma = np.sqrt(self.gp_<%= v.py_varname %>.predict_variances(inputs_gp))
Expand Down Expand Up @@ -75,7 +75,7 @@ class <%= disc.py_classname %>RandomDiscipline(om.ExplicitComponent):

for name in sorted(inputs.keys()):
if name not in random_inputs:
inputs_gp = np.concatenate((inputs_gp, inputs[name]))
inputs_gp = np.concatenate((inputs_gp, inputs[name].flatten()))
inputs_gp = np.atleast_2d(inputs_gp)

<%- disc.output_variables.each do |v| -%>
Expand Down
167 changes: 167 additions & 0 deletions app/lib/whats_opt/templates/egmdo/random_vec_analysis.py.erb
Original file line number Diff line number Diff line change
@@ -0,0 +1,167 @@
import os
import numpy as np
from numpy import nan, inf
import pickle

import openmdao.api as om
<% if @impl.nonlinear_solver.reckless? -%>
from openmdao_extensions.reckless_nonlinear_block_gs import RecklessNonlinearBlockGS
<% else -%>
from openmdao.api import <%= @impl.nonlinear_solver.name %>
<% end -%>
from openmdao.api import NewtonSolver
from openmdao.api import <%= @impl.linear_solver.name %>
<%- @mda.egmdo_random_disciplines.each do |disc| -%>

class <%= disc.py_classname %>RandomVecDiscipline(om.ExplicitComponent):
""" An OpenMDAO base component to encapsulate <%= disc.py_classname %> discipline """
def __init__(self, gp_factory, n_cases=1, **kwargs):
super().__init__(**kwargs)
<%- disc.output_variables.each do |v| -%>
with open(gp_factory.gp_filename("<%= disc.snake_modulename %>", "<%= v.py_varname %>"), 'rb') as f:
self.gp_<%= v.py_varname %> = pickle.load(f)
<%- end %>
self.n_cases = n_cases

<%- unless disc.variables.empty? -%>
def setup(self):
<%- disc.input_variables.numeric.each do |var| %>
<%- if @mda.is_egmdo_random_variable?(var) -%>
<%= var.py_varname %> = np.ones((self.n_cases, <%= var.dim %>))
<%- else -%>
<%= var.py_varname %> = <%= var.init_py_value %>
<%- end -%>
self.add_input('<%= var.py_varname %>', val=<%= var.py_varname %>, desc='<%= var.escaped_desc %>'<%= @impl.use_units && !var.units.blank? ? ", units='#{var.units}'":"" %>)<%- end %>
<%- disc.output_variables.each do |var| %>
<%- if @mda.is_egmdo_random_variable?(var) -%>
_xi_<%= var.py_varname %> = np.zeros(self.n_cases)
<%- else -%>
_xi_<%= var.py_varname %> = 0.0
<%- end -%>
self.add_input('_xi_<%= var.py_varname %>', val=_xi_<%= var.py_varname %>)<%- end %>
<%- disc.output_variables.numeric.each do |var| -%>
<%= var.py_varname %> = np.ones((self.n_cases, <%= var.dim %>))
<%- if var.scaling.blank? -%>
self.add_output('<%= var.py_varname %>', val=<%= var.py_varname %>, desc='<%= var.escaped_desc %>'<%= @impl.use_units && !var.units.blank? ? ", units='#{var.units}'":"" %>)
<%- else -%>
self.add_output('<%= var.py_varname %>', val=<%= var.py_varname %>, desc='<%= var.escaped_desc %>',
ref=<%= var.scaling_ref_py_value %>, ref0=<%= var.scaling_ref0_py_value %>, res_ref=<%= var.scaling_res_ref_py_value %><%= @impl.use_units && !var.units.blank? ? ", units='#{var.units}'":"" %>)
<%- end -%>
<%- end -%>

def compute(self, inputs, outputs):
inputs_gp = np.zeros((self.n_cases, <%= disc.input_variables.map(&:dim).inject(0, :+) %>))
random_inputs = [<%= disc.output_variables.map{|v| "'_xi_#{v.py_varname}'"}.join(', ') %>]
design_vars = [<%= @mda.design_variables.map{|v| "'#{v.name}'"}.join(', ') %>]
coupling_vars = [<%= @mda.coupling_variables.map{|v| "'#{v.name}'"}.join(', ') %>]
i = 0
for name in sorted(inputs.keys()):
if name not in random_inputs:
if name in design_vars:
# design variables should be reshaped
dim = inputs[name].shape[0]
input_temp = np.tile(inputs[name], (self.n_cases, 1))
inputs_gp[:, i:i+dim] = input_temp
i = i + dim
elif name in coupling_vars:
# coupling variables
dim = inputs[name].shape[1]
inputs_gp[:, i:i+dim] = inputs[name]
i = i + dim

inputs_gp = np.atleast_2d(inputs_gp)

<% disc.output_variables.each do |v|%>
sigma = np.sqrt(self.gp_<%= v.py_varname %>.predict_variances(inputs_gp))
mean = self.gp_<%= v.py_varname %>.predict_values(inputs_gp)
outputs['<%= v.py_varname %>'] = mean[:, 0] + inputs['_xi_<%= v.py_varname %>'] * sigma[:, 0]
<%- end -%>
<%- end -%>
<%- end -%>
<%- @mda.plain_disciplines.each do |disc| -%>
<% if !@mda.egmdo_random_disciplines.include?(disc) %>

class <%= disc.py_classname %>RandomVecDiscipline(om.ExplicitComponent):
""" An OpenMDAO base component to encapsulate <%= disc.py_classname %> discipline """
def __init__(self, disc, n_cases=1, **kwargs):
super().__init__(**kwargs)
self.disc = disc
self.n_cases = n_cases

<%- unless disc.variables.empty? -%>
def setup(self):
<%- disc.input_variables.numeric.each do |var| %>
<%- if @mda.is_egmdo_random_variable?(var) -%>
<%= var.py_varname %> = np.ones((self.n_cases, <%= var.dim %>))
<%- else -%>
<%= var.py_varname %> = <%= var.init_py_value %>
<%- end -%>
self.add_input('<%= var.py_varname %>', val=<%= var.py_varname %>, desc='<%= var.escaped_desc %>'<%= @impl.use_units && !var.units.blank? ? ", units='#{var.units}'":"" %>)<%- end %>
<%- disc.output_variables.numeric.each do |var| -%>
<%= var.py_varname %> = np.ones((self.n_cases, <%= var.dim %>))
<%- if var.scaling.blank? -%>
self.add_output('<%= var.py_varname %>', val=<%= var.py_varname %>, desc='<%= var.escaped_desc %>'<%= @impl.use_units && !var.units.blank? ? ", units='#{var.units}'":"" %>)
<%- else -%>
self.add_output('<%= var.py_varname %>', val=<%= var.py_varname %>, desc='<%= var.escaped_desc %>',
ref=<%= var.scaling_ref_py_value %>, ref0=<%= var.scaling_ref0_py_value %>, res_ref=<%= var.scaling_res_ref_py_value %><%= @impl.use_units && !var.units.blank? ? ", units='#{var.units}'":"" %>)
<%- end -%>
<%- end -%>

def compute(self, inputs, outputs):
self.disc.compute(inputs, outputs)

<%- end -%>
<%- end -%>
<%- end -%>


class <%= @mda.py_classname %>RandomVecAnalysis(<%= @impl.parallel_group ? "om.ParallelGroup" : "om.Group" %>):
""" An OpenMDAO base component to encapsulate <%= @mda.py_classname %> random MDA """
def __init__(self, discipline_factory, gp_factory, n_cases=1, **kwargs):
super(). __init__(**kwargs)
self.disc_factory = discipline_factory
self.gp_factory = gp_factory
self.n_cases = n_cases

# self.nonlinear_solver = NewtonSolver(solve_subsystems=False)
self.nonlinear_solver = NonlinearBlockGS()

<% unless @impl.nonlinear_solver.runonce? -%>
self.nonlinear_solver.options['atol'] = <%= @impl.nonlinear_solver.atol %>
self.nonlinear_solver.options['rtol'] = <%= @impl.nonlinear_solver.rtol %>
self.nonlinear_solver.options['err_on_non_converge'] = <%= @impl.to_code(:nonlinear_solver, :err_on_non_converge) %>
self.nonlinear_solver.options['iprint'] = <%= @impl.nonlinear_solver.iprint %>
self.nonlinear_solver.options['maxiter'] = <%= @impl.nonlinear_solver.maxiter %>
<% end -%>

self.linear_solver = <%= @impl.linear_solver.name %>()
self.linear_solver.options['atol'] = <%= @impl.linear_solver.atol %>
self.linear_solver.options['rtol'] = <%= @impl.linear_solver.rtol %>
self.linear_solver.options['err_on_non_converge'] = <%= @impl.to_code(:linear_solver, :err_on_non_converge) %>
self.linear_solver.options['iprint'] = <%= @impl.linear_solver.iprint %>
self.linear_solver.options['maxiter'] = <%= @impl.linear_solver.maxiter %>

def setup(self):
<%- @mda.input_variables.each do |dv| -%>
self.set_input_defaults('<%= dv.name %>', val=<%= dv.init_py_value %><%= @impl.use_units && !dv.units.blank? ? ", units='#{dv.units}'":"" %>)
<%- end -%>
<%- @mda.egmdo_random_variables.each do |v| -%>
self.set_input_defaults('_xi_<%= v.name %>', val=np.zeros(self.n_cases))
<%- end -%>
<%- @mda.plain_disciplines.each do |disc| -%>
name = '<%= disc.py_classname %>RandomVec'
<%- if disc.openmdao_impl&.egmdo_surrogate -%>
disc = <%= disc.py_classname %>RandomVecDiscipline(self.gp_factory, self.n_cases)
<%- else -%>
true_disc = self.disc_factory.create_<%= disc.snake_modulename %>()
disc = <%= disc.py_classname %>RandomVecDiscipline(true_disc, self.n_cases)
<%- end -%>
self.add_subsystem(name, disc, promotes=['*'])
<%- end -%>

Loading

0 comments on commit 2d2d251

Please sign in to comment.