Balance of Plant system

The startup transient and load follow input file have the same structure and only differ by their initial condition and control system governing the transient. The components are the same.

Fluid flow

The FluidProperties blocks are written in the input file in order to declare the properties of the gas in each loop.

[FluidProperties]
  [he]
    type = IdealGasFluidProperties
    molar_mass = 4e-3
    gamma = 1.67
    k = 0.2556
    mu = 3.22639e-5
  []

  [air]
    type = IdealGasFluidProperties
    molar_mass = 29e-3
    gamma = 1.4
    k = 0.025 # W/(m.K)
    mu = 1.8e-5 # Pa.s
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

The friction factor and heat transfer coefficient in each loop are defined using the Closure1PhaseTHM closures and the default options:

[Closures]
  [thm]
    type = Closures1PhaseTHM
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Solid Materials

The solid materials of the many heat structures are declared using a SolidProperties block. It gives the properties of three materials: the H-451 graphite and the fuel, used in the core, and the steel, used in the heat exchanger and in the recuperator. Each solid material is defined using the ThermalFunctionSolidProperties type. Constant densities, thermal conductivities and specific heats are provided.

[SolidProperties]
  [graphite]
    type = ThermalFunctionSolidProperties
    rho = 2160 # kg/m3
    k = 40 # W/(m.K)
    cp = 2100 # J/(kg.K)         approximate mean specific heat of graphite between 800 K (coolant) and 1400 K (fuel)
  []
  [fuel]
    type = ThermalFunctionSolidProperties
    rho = 10970 # kg/m3
    k = 5 # W/(m.K)
    cp = 300 # J/(kg.K)
  []
  [steel]
    type = ThermalFunctionSolidProperties
    rho = 8050 # kg/m3
    k = 45 # W/(m.K)
    cp = 466 # J/(kg.K)
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Components

Geometry

The geometry is defined using the Components system. Each component is defined using a position, orientation and length for 1D components. A width is also defined for 2D components (heat structures). Channels are defined using the FlowChannel1Phase components and heat structures are defined using the HeatStructureCylindrical or HeatStructurePlate depending on the geometry.

Initial conditions

A common initial temperature is chosen for the whole system. It is done to avoid irregular effects in the first few seconds of the simulation due to a high temperature difference between the primary and secondary sides in the heat exchanger (hx). The only exceptions are the sec_pipe1 and hot_leg of the secondary loop, which are respectively connected to the inlet and outlet of the cycle. Their initial temperature is 300 K.

In the primary loop, a 90 bar pressure is defined and is equal to the steady state pressure in that loop. The secondary loop is initialized to the atmospheric pressure.

The velocities are initialized at a very small non-zero value.

Channels

All the channels are defined using FlowChannel1Phase blocks. Depending on the loop, helium or air is declared as the fluid. The initial conditions are those chosen for each loop.

Primary loop components

commentnote

To make it easier to see the space variations of the pressure, mass flow rate, temperature and power through the primary loop using the "PlotOverLine" filter in Paraview, all its components, except the pri_pipe6, are on the same axis and in the same direction. The pri_pipe6 has the opposite direction and is used to close the cycle.

Core

The core geometry is simplified to be able to use a 1D-2D representation. Only one coolant channel is used representing all the coolant channels in the core. It is coupled to a single heat structure representing the moderator and fuel as shown in Figure 1.

This structure is defined using a HeatStructureCylindrical component. The radius of the coolant channel is defined in the design description, and the graphite and fuel outer radii are defined to preserve to their total volume. The average quantity of graphite and of fuel per coolant channel is then computed and used to define the simplified geometry of the heat structure: a cylindrical heat structure, with a coolant cylinder in the middle of a hollow graphite cylinder, itself in a hollow fuel cylinder. This structure is replicated for each coolant channels using the num_rods parameter.

The geometrical parameters of the coolant itself are adapted to represent all the coolant channels:

  • the channel section A is equal to the sum of all the sections of the real coolant channels, because it conducts the same mass flow rate than all the real coolant channels,

  • the channel heating perimeter P_hf is equal to the sum of the perimeters of all the real channels,

  • the hydraulic diameter D_h is the one of a real coolant channel because it is linked with the frictions and pressure drop.

Figure 1: Method used to simplify the core

[Components]
  [core]

    [coolant_channel]
      type = FlowChannel1Phase
      position = '${pri_x_core} ${pri_y_core} 0.'
      orientation = '1 0 0'
      length = ${core_length_channel}
      n_elems = ${core_channel_n_elems}
      A = '${fparse pi * core_nb_coolant_tot * core_radius_coolant * core_radius_coolant}'
      D_h = '${fparse 2 * core_radius_coolant}'

      fp = he
      initial_p = ${pri_press}
      initial_vel = ${vel_ini_pri}
    []

    [hs]
      type = HeatStructureCylindrical
      position = '${pri_x_core} ${pri_y_core} 0.'
      orientation = '1 0 0'
      length = ${core_length_channel}
      n_elems = ${core_channel_n_elems}

      inner_radius = ${core_radius_coolant}
      num_rods = ${core_nb_coolant_tot}
      initial_T = 400

      names = 'graphite_layer fuel_layer'
      widths = '${core_radius_equiv_mod} ${core_radius_equiv_complete_hs}'
      solid_properties = 'graphite fuel'
      solid_properties_T_ref = '0 0' # These materials are independent of temp.
      n_part_elems = '3 3'
      offset_mesh_by_inner_radius = true
    []

    [core_heating]
      type = HeatSourceFromTotalPower
      hs = core/hs
      regions = fuel_layer
      power = total_power
    []

    [core_ht]
      type = HeatTransferFromHeatStructure1Phase
      flow_channel = core/coolant_channel
      hs = core/hs
      hs_side = inner
      P_hf = '${fparse pi * 2 * core_radius_coolant * core_nb_coolant_tot}'
    []
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

A power of 15 MWth is declared using a power block of TotalPower type.

[Components]
  [total_power]
    type = TotalPower
    power = ${tot_power}
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Pressurizer

A pressurizer is added to maintain a 90 bar pressure during the whole simulation, using a InletStagnationPressureTemperature1Phase component. It is located between the core and the heat exchanger, in order to have a pressure as close as possible to 90 bar in the core.

[Components]
  [pressu]
    [pri_jct_2_3_prz]
      type = VolumeJunction1Phase
      connections = 'pri_pipe2:out pri_pipe3:in pressu/pipe_prz:in'
      position = '${pri_x_pipe3} ${pri_y_pipe3} 0.'
      volume = 1e-3
      initial_p = ${pri_press}
      initial_vel_x = ${vel_ini_pri}
      use_scalar_variables = false
    []
    [pipe_prz]
      type = FlowChannel1Phase
      position = '${pri_x_pipe3} ${pri_y_pipe3} 0.'
      orientation = '0 1 0'
      length = ${PRI_L_prz}
      n_elems = ${pri_pipe_prz_n_elems}
      A = '${pri_pipes_area}'
      D_h = '${pri_pipes_D_h}'

      fp = he
      initial_p = ${pri_press}
      initial_vel = ${vel_ini_pri}
    []
    [prz]
      type = InletStagnationPressureTemperature1Phase
      p0 = ${pri_press}
      T0 = ${T_ini}
      input = 'pressu/pipe_prz:out'
    []
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Heat exchanger

The heat exchanger is modeled as a 2 meters-high structure, composed of couples of primary and secondary channels (see Figure 2 ) which are replicated 20,000 times. Their diameters are respectively 3 mm and 5 mm. They are defined using FlowChannel1Phase components. They are separated by a 1 mm steel wall defined using a HeatStructureCylindrical component, and the flow directions are opposite. The heat transfers are declared with HeatTransferFromHeatStructure1Phase. The wall is used as heat structure and each channel as the heated flow channel.

The same principle as for the core coolant channel is used to define the A section, P_hf heating perimeter and D_h hydraulic diameter parameters. The replication of the couple is done using the num_rods parameter of the wall block.

Figure 2: Couple of primary and secondary flow channels in the heat exchanger

[Components]
  [hx]
    [pri]
      type = FlowChannel1Phase
      position = '${pri_x_hx} ${pri_y_hx} 0.'
      orientation = '1 0 0'
      length = ${hx_length}
      n_elems = ${hx_n_elems_axial}
      A = '${fparse pi * hx_nb_channels * hx_dia_pri * hx_dia_pri * 0.25}'
      D_h = ${hx_dia_pri}

      fp = he
      initial_p = ${pri_press}
      initial_vel = ${vel_ini_pri}
    []

    [ht_pri]
      type = HeatTransferFromHeatStructure1Phase
      hs = hx/wall
      hs_side = inner
      flow_channel = hx/pri
      P_hf = '${fparse pi * hx_nb_channels * hx_dia_pri}'
    []

    [wall]
      type = HeatStructureCylindrical
      length = ${hx_length}
      n_elems = ${hx_n_elems_axial}
      n_part_elems = 1
      names = 'hx_wall'
      orientation = '1 0 0'
      position = '${pri_x_hx} ${pri_y_hx} 0.'
      widths = '${hx_wall_thickness}'
      solid_properties = 'steel'
      solid_properties_T_ref = '0' # This material is independent of temp.
      inner_radius = '${fparse hx_dia_pri / 2}'
      num_rods = ${hx_nb_channels}
    []

    [ht_sec]
      type = HeatTransferFromHeatStructure1Phase
      hs = hx/wall
      hs_side = outer
      flow_channel = hx/sec
      P_hf = '${fparse pi * hx_nb_channels * hx_dia_sec}'
    []

    [sec]
      type = FlowChannel1Phase
      position = '${sec_x_hx} ${sec_y_hx} 0.'
      orientation = '-1 0 0'
      length = ${hx_length}
      n_elems = ${hx_n_elems_axial}
      A = '${fparse pi * hx_nb_channels * hx_dia_sec * hx_dia_sec * 0.25}'
      D_h = '${hx_dia_sec}'

      fp = air
      initial_p = ${p_sec}
      initial_vel = ${vel_ini_sec}
    []
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Pump

The pressure drop is compensated by a pump declared as a ShaftConnectedPump1Phase component. It is connected to a motor, declared as ShaftConnectedMotor, and a shaft, declared as Shaft, which provide its hydraulic torque. The rated values of the pump (volumetric flow rate, head, density, torque, shaft speed) are defined as close as possible to the operating values (Table 1).

Table 1: Pump - comparison between the rated and steady state values

ParameterRated valueSteady state value
Density (kg/m3)54.86
Volumetric flow rate (m3/s)21.92
Pump head (m)350293
Shaft speed (rad/s)54.5
Torque (N.m)5042

The performance curves (the Bingham head and torque functions) are defined using a head_fcn and a torque_fcn functions in the Functions block, which imports data from csv files. The motor torque is set to match the nominal mass flow rate in the loop.

[Components]
  [circ]
    [pump]
      type = ShaftConnectedPump1Phase
      inlet = 'pri_pipe4:out'
      outlet = 'pri_pipe5:in'
      position = '${pri_x_pipe5} ${pri_y_pipe5} 0.'
      A_ref = ${pump_area}
      scaling_factor_rhoEV = 1e-5
      volume = ${pump_volume}
      inertia_coeff = ${pump_inertia_coeff}
      inertia_const = ${pump_inertia_const}
      omega_rated = ${pump_omega_rated}
      speed_cr_I = ${pump_speed_cr_I}
      speed_cr_fr = ${pump_speed_cr_fr}
      torque_rated = ${pump_torque_rated}
      volumetric_rated = ${pump_volumetric_rated}
      head_rated = ${pump_head_rated}
      tau_fr_coeff = ${pump_tau_fr_coeff}
      tau_fr_const = ${pump_tau_fr_const}
      head = head_fcn
      torque_hydraulic = torque_fcn
      density_rated = ${pump_density_rated}
      initial_p = ${pri_press}
      initial_vel_x = ${vel_ini_pri}
      use_scalar_variables = false
    []

    [motor]
      type = ShaftConnectedMotor
      inertia = ${pri_motor_inertia}
      torque = ${pri_motor_torque}
    []

    [shaft]
      type = Shaft
      connected_components = 'circ/motor circ/pump'
      initial_speed = ${shaft_initial_speed}
    []
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Secondary loop components

The secondary loop is a recuperated open-air Brayton cycle and is similar to the one described in the Thermal-Hydraulics Module Modeling Guide.

This input file has been adapted to achieve the nominal operating condition. In particular, the rated values for the turbine and the compressor and the geometrical parameters were adjusted.

The detail of the components is given below.

Boundary Conditions

Stagnation pressure and temperature are defined at the inlet of the Brayton cycle using an InletStagnationPressureTemperature1Phase component with = 300 K and = 1 bar.

[Components]
  [inlet]
    type = InletStagnationPressureTemperature1Phase
    input = 'sec_pipe1:in'
    p0 = ${p_ambient}
    T0 = ${T_ambient}
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

The pressure is prescribed at the outlet using an Outlet1Phase component.

[Components]
  [outlet]
    type = Outlet1Phase
    input = 'hot_leg:out'
    p = ${p_ambient}
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Shaft

The compressor, the turbine, the generator, and the motor are on the same shaft. This component is defined using a Shaft type.

[Components]
  [shaft]
    type = Shaft
    connected_components = 'motor compressor turbine generator'
    initial_speed = ${speed_initial}
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Compressor

The compressor is defined using the ShaftConnectedCompressor1Phase component type. The A_ref section is adapted to admit enough air in the cycle. The rated values (omega_rated, mdot_rated, c0_rated, rho0_rated, see Figure 3) are defined to match as close as possible to the operating values. The results are given, with those of the turbine, in the next table.

These rated parameters are associated to the efficiency and pressure ratio functions of the compressor, defined in the Functions block and using data from csv files.

[Components]
  [compressor]
    type = ShaftConnectedCompressor1Phase
    position = '${sec_x_pipe2} ${sec_y_pipe2} 0'
    inlet = 'sec_pipe1:out'
    outlet = 'sec_pipe2:in'
    A_ref = ${A_ref_comp}
    volume = ${V_comp}

    omega_rated = ${speed_rated}
    mdot_rated = ${rated_mfr}
    c0_rated = ${c0_rated_comp}
    rho0_rated = ${rho0_rated_comp}

    # Determines which compression ratio curve and efficiency curve to use depending on ratio of speed/rated_speed
    speeds = '0.5208 0.6250 0.7292 0.8333 0.9375'
    Rp_functions = 'rp_comp1 rp_comp2 rp_comp3 rp_comp4 rp_comp5'
    eff_functions = 'eff_comp1 eff_comp2 eff_comp3 eff_comp4 eff_comp5'

    min_pressure_ratio = 1.0

    speed_cr_I = 0
    inertia_const = ${I_comp}
    inertia_coeff = '${I_comp} 0 0 0'

    # assume no shaft friction
    speed_cr_fr = 0
    tau_fr_const = 0
    tau_fr_coeff = '0 0 0 0'

    initial_p = ${p_sec}
    initial_vel_x = ${vel_ini_sec}

    use_scalar_variables = false
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Turbine

A turbine is defined using the ShaftConnectedCompressorTurbine component type: a turbine can be considered as an inverted compressor. The treat_as_turbine parameter is defined as “true”.

Once again, the typical section defined by A_ref and the rated values (omega_rated, mdot_rated, c0_rated, rho0_rated, see Figure 3) are adapted to match as close as possible with the operating values. The results are given in Table 2.

Table 2: Compressor and turbine - comparison between the rated and steady state values

CompressorTurbine
Rated valuesSteady state valuesRated valuesSteady state values
Density (kg/m3)1.20.81.41.43
Shaft speed (rad/s)9948932999489329
Sound speed (m/s)340330670676
Mass flow rate (kg/s)2020.52020.5

These rated parameters are associated with the efficiency and pressure ratio functions of the turbine, defined in the Functions block and using data from csv files.

Figure 3: Compressor and turbine rated parameters

The shaft speed is smaller than the rated value. This difference is due to a margin between the rated value and the value imposed by the PID controller (presented in the next part). In a real system, the shaft speed should not be over the rated value, to avoid damage to the components.

[Components]
  [turbine]
    type = ShaftConnectedCompressor1Phase
    position = '${sec_x_pipe6} ${sec_y_pipe6} 0'
    inlet = 'sec_pipe5:out'
    outlet = 'sec_pipe6:in'
    A_ref = ${A_ref_turb}
    volume = ${V_turb}

    # A turbine is treated as an "inverse" compressor, this value determines if component is to be treated as turbine or compressor
    # If treat_as_turbine is omitted, code automatically assumes it is a compressor
    treat_as_turbine = true

    omega_rated = ${speed_rated}
    mdot_rated = ${rated_mfr}
    c0_rated = ${c0_rated_turb}
    rho0_rated = ${rho0_rated_turb}

    # Determines which compression ratio curve and efficiency curve to use depending on ratio of speed/rated_speed
    speeds = '0 0.5208 0.6250 0.7292 0.8333 0.9375'
    Rp_functions = 'rp_turb0 rp_turb1 rp_turb2 rp_turb3 rp_turb4 rp_turb5'
    eff_functions = 'eff_turb1 eff_turb1 eff_turb2 eff_turb3 eff_turb4 eff_turb5'

    min_pressure_ratio = 1.0

    speed_cr_I = 0
    inertia_const = ${I_turb}
    inertia_coeff = '${I_turb} 0 0 0'

    # assume no shaft friction
    speed_cr_fr = 0
    tau_fr_const = 0
    tau_fr_coeff = '0 0 0 0'

    initial_p = ${p_sec}
    initial_vel_x = ${vel_ini_sec}
    use_scalar_variables = false
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Generator

A generator is added on the same shaft as the turbine to produce electricity. The ShaftConnectedMotor component type is once again used. The generator torque is defined as follow:

C is defined to generate 2 MWe at the rated shaft speed:

This parameter is used in the Functions block by the generator_torque_fn. This one is used to define the torque of the generator in the component block. Note that a negative sign is present because energy is extracted.

[Components]
  [generator]
    type = ShaftConnectedMotor
    inertia = ${I_generator}
    torque = generator_torque_fn
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Motor

A motor component of ShaftConnectedMotor type is added and connected to the Shaft. Its behavior is defined by the Control and ControlLogic blocks (see below).

[Components]
  [motor]
    type = ShaftConnectedMotor
    inertia = ${I_motor}
    torque = 0 # controlled
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Recuperator

The recuperator includes two counter-current coolant channels (the hot_leg and cold_leg). A HeatStructureCylindrical component is added to couple these two channels, and two heat transfers are defined using the HeatTransferFromHeatStructure1Phase. Each of them couples a channel with the heat structure.

[Components]
  [cold_leg]
    type = FlowChannel1Phase
    position = '${sec_x_cold} ${sec_y_cold} 0'
    orientation = '0 -1 0'
    length = '${SEC_L_COLD}'
    n_elems = '${sec_n_elems_cold}'
    A = ${SEC_A_COLD}

    fp = air
    initial_p = ${p_sec}
    initial_vel = ${vel_ini_sec}
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)
[Components]
  [hot_leg]
    type = FlowChannel1Phase
    position = '${sec_x_hot} ${sec_y_hot} 0'
    orientation = '0 1 0'
    length = ${SEC_L_HOT}
    n_elems = ${sec_n_elems_hot}
    A = ${SEC_A_HOT}

    fp = air
    initial_T = ${T_ambient}
    initial_p = ${p_ambient}
    initial_vel = ${vel_ini_sec}
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)
[Components]
  [recuperator]
    type = HeatStructureCylindrical
    orientation = '0 -1 0'
    position = '${sec_x_cold} ${sec_y_cold} 0'
    length = '${SEC_L_COLD}'
    widths = ${recuperator_width}
    n_elems = '${sec_n_elems3}'
    n_part_elems = 2
    names = recuperator
    solid_properties = steel
    solid_properties_T_ref = '0' # This material is independent of temp.
    inner_radius = '${fparse SEC_D1 / 2}'
    offset_mesh_by_inner_radius = true
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)
[Components]
  [heat_transfer_cold_leg]
    type = HeatTransferFromHeatStructure1Phase
    flow_channel = cold_leg
    hs = recuperator
    hs_side = OUTER
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)
[Components]
  [heat_transfer_hot_leg]
    type = HeatTransferFromHeatStructure1Phase
    flow_channel = hot_leg
    hs = recuperator
    hs_side = INNER
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Controls

A ControlLogic block is added to control the motor torque.

The initial_motor_PID of type PIDControl defines the behavior of the motor at the beginning. It depends on the set_point block, whose type is a GetFunctionValueControl. It declares the desired shaft speed, which is the rated shaft speed minus 9000 rpm (see the turbine paragraph)

As the shaft speed and heat transfer from the core progressively increase, the turbine generates torque and the need for the motor torque becomes smaller. Thus, it progressively decreases, and it is eventually shut off.

A motor_PID, whose type is SetComponentRealValueControl, decides which value must be attributed to the motor torque. If the motor torque is too small, the motor is shut, using a motor_torque_fn_shutdown function written in the Functions block. Else, the behavior of the motor is still defined by the PIDControl. This condition is defined using the logic block, whose type is ParsedFunctionControl. In this one, a time condition is added: if the elapsed time is too small, the condition on the motor torque value is not applicable. At the beginning of the simulation: the motor torque is small and increasing, but it should not be shut.

[ControlLogic]

  # Sets desired shaft speed to be reached by motor NOTE: SHOULD BE SET LOWER THAN RATED TURBINE RPM
  [set_point]
    type = GetFunctionValueControl
    function = '${fparse speed_rated_rpm - 9000}'
  []

  # PID with gains determined by iterative process NOTE: Gain values are system specific
  [initial_motor_PID]
    type = PIDControl
    set_point = set_point:value
    input = shaft_RPM

    initial_value = ${motor_ini_val}
    K_p = ${motor_K_p}
    K_i = ${motor_K_i}
    K_d = ${motor_K_d}
  []

  # Determines when the PID system should be running and when it should begin the shutdown cycle. If needed: PID output, else: shutdown function
  [logic]
    type = ParsedFunctionControl
    function = 'if(time_simulation < 100, PID, if(motor > ${motor_shutdown_value}, PID, shutdown_fn))'
    symbol_names = 'time_simulation motor turb PID shutdown_fn'
    symbol_values = 'time_simulation sec_motor_torque turbine_torque initial_motor_PID:output sec_motor_torque_fn_shutdown'
  []

  # Takes the output generated in [logic] and applies it to the motor torque
  [motor_PID]
    type = SetComponentRealValueControl
    component = motor
    parameter = torque
    value = logic:value
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)

Startup and load follow transient

The model description above defines mainly the model used for the startup transient. In the case of the load follow transient, some blocks and parameters are disabled, because the system starts from the operating values computed in the startup transient. The following elements are commented:

  • the initial_vel, initial_vel_x, initial_vel_y, initial_vel_z, initial_T, initial_p in the GlobalParams and in the various components where they were defined (pipes, junctions, etc.),

  • the is_tripped_fn, PID_tripped_constant_value, PID_tripped_status_fn, sec_motor_torque_fn_shutdown Functions, which defined in the startup transient the motor behavior,

  • the time_trip_aux and PID_trip_status_aux AuxScalarKernels, which associated the previous functions results to the time_trip and PID_trip_status AuxVariables,

  • the ControlLogic blocks that gave the retroaction on the motor behavior: the set_point, initial_motor_PID, logic and motor_PID,

  • the Controls that are associated: PID_trip_status and time_PID.

The AuxVariables defined to control the motor behavior are not disabled because the input file of the load follow transient needs to import data from the startup transient. To do so, some of the defined objects in the two input file, particularly the AuxVariables, need to be the same.

Two ControlLogic blocks are added to impose the core power during the load follow transient. A power_logic, whose type is ParsedFunctionControl, calls a power_fn function, which defines the power dependance on the time. It is then transferred to a power_applied sub block (of type SetComponentRealValueControl) which indicates to attribute this behavior to the total_power component block.

[ControlLogic]

  # Changes of power
  [power_logic]
    type = ParsedFunctionControl
    function = 'power_fn'
    symbol_names = 'power_fn'
    symbol_values = 'power_fn'
  []
  # Applies heat source to the total_power block
  [power_applied]
    type = SetComponentRealValueControl
    component = total_power
    parameter = power
    value = power_logic:value
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_load_follow_transient.i)

The data are transferred from the startup transient results to the load follow simulation using a supplemental output in the first input file. In the Outputs block, the checkpoint parameter is defined as true. The results saved in this output are then called in a Problem block in the load follow input file.

[Outputs]

  file_base = 'htgr_startup_transient'
  exodus = true
  [csv]
    type = CSV
  []

  [cp]
    type = Checkpoint
    execute_on = 'FINAL'
  []

  [console]
    type = Console
    max_rows = 5
    show = 'core_T_in core_T_out core_p_difference core_m_dot_in core_power_difference
    hx_pri_T_in hx_pri_T_out  hx_pri_p_difference  hx_pri_power_difference
    pump_p_in pump_p_out pump_p_difference

    comp_T_in comp_T_out comp_p_in comp_p_out
    turb_T_in turb_T_out turb_p_in turb_p_out
    generator_power cycle_efficiency shaft_speed
    hx_sec_T_in hx_sec_T_out  hx_sec_p_in hx_sec_p_out hx_sec_p_difference hx_sec_m_dot_in
    hx_sec_m_dot_difference hx_sec_power_flux
    '
  []
[]
(microreactors/gcmr/balance_of_plant/htgr_startup_transient.i)
[Problem]
  restart_file_base = htgr_startup_transient_cp/LATEST
  # Component areas are set using initial conditions
  allow_initial_conditions_with_restart = true
[]
(microreactors/gcmr/balance_of_plant/htgr_load_follow_transient.i)