How to use the openmdao.api.Problem function in openmdao

def test_BladeCurvature():

    data = {}
    data['r'] = np.array([1.5375, 1.84056613137, 1.93905987557, 2.03755361977, 2.14220322299, 2.24069696719, 2.33919071139, 2.90419974, 3.04095863882, 3.13945238302, 5.637500205, 7.04226699698, 8.370800055, 11.8494282928, 13.8375, 14.7182548841, 15.887499795, 18.5537233505, 19.9875, 22.0564071286, 24.087500205, 26.16236509, 28.187499795, 32.2875, 33.5678634821, 36.387500205, 38.5725768593, 40.487499795, 40.6967540173, 40.7964789782, 40.8975, 42.6919644014, 44.5875, 52.787499795, 56.204199945, 58.937499795, 61.67080026, 63.0375])
    data['precurve'] = np.array([0.0, 0.043324361025, 0.0573893371698, 0.0714469497372, 0.0863751069858, 0.100417566593, 0.114452695996, 0.194824077331, 0.214241777459, 0.228217752953, 0.580295739194, 0.776308800624, 0.960411633829, 1.4368012564, 1.7055214864, 1.823777005, 1.98003324362, 2.32762426752, 2.50856911855, 2.76432512112, 3.0113656418, 3.26199245912, 3.50723775206, 4.0150233695, 4.17901272929, 4.55356019347, 4.85962948702, 5.14086873143, 5.17214747287, 5.18708601127, 5.20223968442, 5.47491847385, 5.77007321175, 7.12818875977, 7.7314427824, 8.22913789456, 8.73985955154, 9.0])
    data['presweep'] = np.array([0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0])
    data['precone'] = 2.5
    npts = np.size(data['r'])

    prob = Problem()
    prob.root = Group()
    prob.root.add('comp', BladeCurvature(npts), promotes=['*'])
    prob = init_IndepVar_add(prob, data)
    prob.setup(check=False)
    prob = init_IndepVar_set(prob, data)

    check_gradient_unit_test(prob, tol=5e-5)

except ImportError:
    print("PETSc is not installed. Skipping beam optimization.")
    PETScVector = None


if __name__ == '__main__' and PETScVector is not None:
    E = 1.
    L = 1.
    b = 0.1
    volume = 0.01

    num_cp = 5
    num_elements = 50
    num_load_cases = 2

    prob = om.Problem(model=MultipointBeamGroup(E=E, L=L, b=b, volume=volume,
                                                num_elements=num_elements, num_cp=num_cp,
                                                num_load_cases=num_load_cases))

    prob.driver = om.ScipyOptimizeDriver()
    prob.driver.options['optimizer'] = 'SLSQP'
    prob.driver.options['tol'] = 1e-9
    prob.driver.options['disp'] = True

    prob.setup()

    prob.run_driver()

    h = prob['interp.h']
    expected = np.array([ 0.14122705,  0.14130706,  0.14154096,  0.1419107,   0.14238706,  0.14293095,
                          0.14349514,  0.14402636,  0.1444677,   0.14476123,  0.14485062,  0.14468388,
                          0.14421589,  0.1434107,   0.14224356,  0.14070252,  0.13878952,  0.13652104,

def test_implicit_ext_solve(self):
        top = Problem()
        root = top.root = Group()
        root.add('p1', IndepVarComp('x', 0.5))
        root.add('comp', SimpleImplicitComp())
        root.add('comp2', ExecComp('zz = 2.0*z'))

        root.connect('p1.x', 'comp.x')
        root.connect('comp.z', 'comp2.z')

        root.ln_solver = ScipyGMRES()
        root.nl_solver = Newton()
        top.setup(check=False)

        top.run()
        assert_rel_error(self, top['comp.z'], 2.666667, 1e-5)

def test_two_simple(self):
        size = 3
        group = Group()
        group.add('P', IndepVarComp('x', numpy.ones(size)))
        group.add('C1', DistribExecComp(['y=2.0*x'], arr_size=size,
                                        x=numpy.zeros(size),
                                        y=numpy.zeros(size)))
        group.add('C2', ExecComp(['z=3.0*y'],
                                 y=numpy.zeros(size),
                                 z=numpy.zeros(size)))

        prob = Problem(impl=impl)
        prob.root = group
        prob.root.ln_solver = LinearGaussSeidel()
        prob.root.connect('P.x', 'C1.x')
        prob.root.connect('C1.y', 'C2.y')

        prob.setup(check=False)
        prob.run()

        J = prob.calc_gradient(['P.x'], ['C2.z'], mode='fwd', return_format='dict')
        assert_rel_error(self, J['C2.z']['P.x'], numpy.eye(size)*6.0, 1e-6)

        J = prob.calc_gradient(['P.x'], ['C2.z'], mode='rev', return_format='dict')
        assert_rel_error(self, J['C2.z']['P.x'], numpy.eye(size)*6.0, 1e-6)

def execute(self):

		# --- Import Modules
		import numpy as np
		import os
		from openmdao.api import IndepVarComp, Component, Group, Problem, Brent, ScipyGMRES, ScipyOptimizer, DumpRecorder
		from rotorse.rotor_aeropower import RotorAeroPower
		from rotorse.rotor_geometry import RotorGeometry, NREL5MW, DTU10MW, NINPUT
		from rotorse import RPM2RS, RS2RPM, TURBULENCE_CLASS, DRIVETRAIN_TYPE

		rotor = Problem()
		myref = DTU10MW()

		npts_coarse_power_curve = 20 # (Int): number of points to evaluate aero analysis at
		npts_spline_power_curve = 200  # (Int): number of points to use in fitting spline to power curve

		rotor.root = RotorAeroPower(myref, npts_coarse_power_curve, npts_spline_power_curve)
		# ---

		# --- Optimizer
		rotor.driver = ScipyOptimizer()
		rotor.driver.options['optimizer'] = 'SLSQP'
		rotor.driver.options['tol'] = 1.0e-8
		# ---
		# --- Objective
		# AEP0 = 47147596.2911617
		AEP0 = 48113504.25433461

def main():  # pragma: no cover
    p = om.Problem()
    p.model = LaunchVehicleLinearTangentODE(num_nodes=1)
    p.setup()

    p.run_model()

    var_map = {'rho': r'$\rho$',
               'x': r'$x$',
               'y': r'$y$',
               'vx': r'$v_x$',
               'vy': r'$v_y$',
               'm': r'$m$',
               'mdot': r'$\dot{m}$',
               'Isp': r'$I_{sp}$',
               'thrust': r'$F_{T}$',
               'ydot': r'$\dot{y}$',
               'xdot': r'$\dot{x}$',

'strength_factor_for_upper_skin' : 1.0, # the yield stress is multiplied by this factor for the upper skin

                    'wing_weight_ratio' : 1.25,
                    'exact_failure_constraint' : False, # if false, use KS function

                    'struct_weight_relief' : True,
                    'distributed_fuel_weight' : True,

                    'fuel_density' : 803.,      # [kg/m^3] fuel density (only needed if the fuel-in-wing volume constraint is used)
                    'Wf_reserve' :15000.,       # [kg] reserve fuel mass
                    }

        surfaces = [surf_dict]

        # Create the problem and assign the model group
        prob = Problem()

        # Add problem information as an independent variables component
        indep_var_comp = IndepVarComp()
        indep_var_comp.add_output('v', val=np.array([.85 * 295.07, .64 * 340.294]), units='m/s')
        indep_var_comp.add_output('alpha', val=0., units='deg')
        indep_var_comp.add_output('alpha_maneuver', val=0., units='deg')
        indep_var_comp.add_output('Mach_number', val=np.array([0.85, 0.64]))
        indep_var_comp.add_output('re',val=np.array([0.348*295.07*.85*1./(1.43*1e-5), \
                                  1.225*340.294*.64*1./(1.81206*1e-5)]),  units='1/m')
        indep_var_comp.add_output('rho', val=np.array([0.348, 1.225]), units='kg/m**3')
        indep_var_comp.add_output('CT', val=0.53/3600, units='1/s')
        indep_var_comp.add_output('R', val=14.307e6, units='m')
        indep_var_comp.add_output('W0', val=148000 + surf_dict['Wf_reserve'],  units='kg')
        indep_var_comp.add_output('speed_of_sound', val= np.array([295.07, 340.294]), units='m/s')
        indep_var_comp.add_output('load_factor', val=np.array([1., 2.5]))
        indep_var_comp.add_output('empty_cg', val=np.zeros((3)), units='m')

def execute(self):

		import numpy as np
		import os
		from openmdao.api import IndepVarComp, Component, Group, Problem, Brent, ScipyGMRES
		from rotorse.rotor_aeropower import RotorAeroPower
		from rotorse.rotor_geometry import RotorGeometry, NREL5MW, DTU10MW, NINPUT
		from rotorse import RPM2RS, RS2RPM, TURBULENCE_CLASS, DRIVETRAIN_TYPE, TURBINE_CLASS



		myref = DTU10MW()
		rotor = Problem()

		npts_coarse_power_curve = 20 # (Int): number of points to evaluate aero analysis at
		npts_spline_power_curve = 200  # (Int): number of points to use in fitting spline to power curve

		rotor.root = RotorAeroPower(myref, npts_coarse_power_curve, npts_spline_power_curve)
		rotor.setup()

		# === blade grid ===
		rotor['hubFraction'] = 0.023785  # (Float): hub location as fraction of radius
		rotor['bladeLength'] = 96.7  # (Float, m): blade length (if not precurved or swept) otherwise length of blade before curvature
		rotor['precone'] = 4.  # (Float, deg): precone angle
		rotor['tilt'] = 6.0  # (Float, deg): shaft tilt
		rotor['yaw'] = 0.0  # (Float, deg): yaw error
		rotor['nBlades'] = 3  # (Int): number of blades
		rotor['r_max_chord'] = 0.2366  # (Float): location of max chord on unit radius

surf_dict = default_dict.copy()

    mesh = meshes[surf_name]

    surf_dict['name'] = surf_name
    surf_dict['mesh'] = mesh

    if surf_name == 'Tail_H':
        # Set one twist variable on tail to control incidence for trim
        surf_dict['twist_cp'] = np.zeros(1)

    if surf_name not in dont_include:
        surfaces.append(surf_dict)

# Create the problem and the model group
prob = Problem()

# Define flight variables as independent variables of the model
indep_var_comp = IndepVarComp()
indep_var_comp.add_output('v', val=248.136, units='m/s')
indep_var_comp.add_output('alpha', val=5., units='deg')
indep_var_comp.add_output('beta', val=0., units='deg')
indep_var_comp.add_output('Mach_number', val=0.84)
indep_var_comp.add_output('re', val=1.e6, units='1/m')
indep_var_comp.add_output('rho', val=0.38, units='kg/m**3')
indep_var_comp.add_output('cg', val=np.array([33.68, 0.0, 4.52]), units='m')
indep_var_comp.add_output('S_ref_total', val=383.7, units='m**2')
indep_var_comp.add_output('omega', val=np.array([0.0, 0.0, 0.0]), units='deg/s')

prob.model.add_subsystem('prob_vars',
    indep_var_comp,
    promotes=['*'])

pred_cok = []
    inputs = np.linspace(0, 1, ngrid)
    for x in inputs:
        pbcok['mm.x'] = x
        pbcok.run()
        pred_cok.append(pbcok['mm.f_x'])

    pbcok_mu    = np.array([float(p[0]) for p in pred_cok])
    pbcok_sigma = np.array([float(p[1]) for p in pred_cok])

    ## "Co-kriging" with 1 level of fidelity a.k.a. kriging
    surrogate = MultiFiCoKrigingSurrogate()
    ## Kriging from openmdao
    #surrogate = KrigingSurrogate()

    pbk = Problem(Simulation(surrogate, nfi=1))
    pbk.setup()

    pbk['mm.train:x'] = np.array(doe_e).reshape(len(doe_e),1)
    pbk['mm.train:f_x'] = model_hifi(pbk['mm.train:x'])

    pbk.run() # train

    ngrid = 100
    pred_k = []
    inputs = np.linspace(0, 1, ngrid)
    for x in inputs:
        pbk['mm.x'] = x
        pbk.run()
        pred_k.append(pbk['mm.f_x'])

    pbk_mu    = np.array([float(p[0]) for p in pred_k])