AeroCable.int

/***
 -------------------------------------------------------------------------------
 |                                                                             |
 | NASA Glenn Research Center                                                  |
 | 21000 Brookpark Rd                                                          |
 | Cleveland, OH 44135 	                                                       |
 |                                                                             |
 | File Name:     AeroCable.int                                                |
 | Author(s):     George Thomas, Jeffrey Csank, David Sadey, Tom Lavelle       |
 | Date(s):       April 2018                                                   |
 |                                                                             |
 -------------------------------------------------------------------------------
***/

#ifndef __AEROCABLE__
#define __AEROCABLE__

#include "PSL_math.fnc"

class AeroCable extends ElectricElement {

  //----------------------------
  // ****** DOCUMENTATION ******
  //----------------------------  

  // title = "";

  description = "The " + isA() + " is a cable model that uses lookup table data based on real 
world cable metrics (resistance, inductance, and mass per unit length as well as ampacity), 
as a function of wire gauge. Sizing logic is included and used during on-design runs, that 
selects the lightest possible cable assembly (bundle of n number of parallel cables) that 
meets the ampacity requirements (on-design current).";

  usageNotes = isA() +
  "

- The input design parameters that must be specified by the modeler when creating an instance 
of this " + isA() + " component are cable length (length), maximum number of parallel wires 
(n_parallel) and voltage drop tolerance, expressed as a percentage of the cable input voltage 
(drop_tolerance) When this component is run ONDESIGN, it will select a wire gauge and number 
of parallel wires that comprises a cable assembly rated above the cable current seen during 
this ONDESIGN run. Outputs generated by ONDESIGN runs include the cable resistance and 
inductance (R and L), the cable mass (Mass), and mass specific ampacity (SpecificAmperage).

- This component supports DC and AC (single or multiphase) power (see ElectricPowerType). 
At instantiation, this component's power type is unspecified. NPSS Power System Components 
need to have a power type in order to correctly run power calculations. For this component, 
the power type is typically set by calling the propagateEpsSolverListAndPowerTypes() method at the top 
level assembly. This method is defined in ElectricPort.prt.

- The data used in this design logic was obtained from lookup tables constructed using 
realistic cable data. Note that this data is based on terrestrial medium voltage power cable, 
and not flight-weight, medium voltage aircraft cable. This cable does not yet exist 
(as of 2020), and so be aware that the cable designs produced by this sizing logic may 
overpredict mass.

- This component does not contain a node, and so starts each iteration with its port voltages 
unknown. The concept in NPSS Power System Library is that all components that contain nodes 
all have a prePass() method specified that pushes node voltage information, via their ports
to all components that do not contain nodes and so do not know their voltages. The intent is 
that after all component prePass methods are automatically called, all components in the 
system will then have port voltages specified. Then during the normal iteration, components 
that don't contain nodes such as this one, will execute first. These components typically 
transmit power from one node to another and do not transform it from one kind to another 
(e.g. cables and breaker). These execute first in the solver sequence, and inside their 
calculate() methods, they calculate their currents, which are then passed via the ports, to 
all other components in the circuit. So after these nodeless transmission components 
execute, all electric port voltages and currents should be specified, and so all other 
components are free to run.

- This component, like other power system components in the NPSS Power System Library, can 
optionally include thermal models. An optional thermal model is enabled by setting 
switchThermPort to TRUE, and plugging an EThermalMass subelement into the S_eThermMass socket. 
Doing these will add a temperature state (existing within EThermalMass) and a thermal port to 
the model. The thermal port is intended to connect this component to a second component that 
represents the mechanism by which heat is extracted from this component. This second component 
could represent a heat exchanger, cold plate, or just model heat transfer from the first 
component to the surrounding environment. For more information, see EThermalMass. 
 
- Solver independents and dependents are only used in the " + isA() + " component during 
ONDESIGN. In ONDESIGN, the cable's impedance (resistance and inductance R and L) is not known. 
Because of this, the cable will run with only its input and output voltages specified via its 
input and output ports (see NPSS Power Systems Library - solverSequence). With no other 
parameters, the cable cannot know the design current (ampacity) it must be sized to. For this 
reason, the cable includes solver variables for ONDESIGN. The cable's real current (and also
imaginary current if the cable carries AC) is set up as an independent. A dependent is used
to ensure that the voltage seen across the cable's ports, equals the voltage drop according
to the guessed current and designed impedance (dV = dVCalc = I * Z). Imaginary dependents
are also included if the cable power type is AC (AC1 or AC3).

- The " + isA() + " OFFDESIGN logic does not need solver variables because during OFFDESIGN, 
the cable knows its complex impedance (Z = R + j*2*PI*frequency*L), and knows its voltage as 
specified at its input and output ports. Current is calculated from I = dV/Z, and this 
current information is updated in the cable's input and output ports.";
  
  //------------------------------
  // ****** SETUP VARIABLES ******
  //------------------------------
  
  real R {
    value = .01;  IOstatus = "input";  units = "ohm";
  }
  
  real L {
    value = .0000796;  IOstatus = "input";  units = "H"; // [henrys]
  }
  
  real length {
    value = 1;  IOstatus = "input";  units = "m"; // [meters]
  }

  // will add losses in future releases.
  Loss_r {
    value = 0; IOstatus = "output"; units = "W"; // [watts]
    description = "Cable loss.";
  }

  Loss_j {
    value = 0; IOstatus = "output"; units = "W"; // [watts]
    description = "Cable loss.";
  }

  Q_heat {
    value = 0; IOstatus = "output"; units = "Btu/sec"; // [BTUs / sec]
    description = "Power dissipation at current time.";
  }

  Mass {
    value = 1; IOstatus = "output"; units = "kg"; // [kilograms]
    description = "Mass output of the cable.";
  }
  
  real SpecificAmperage {
    value = 170; IOstatus = "input"; units = "none"; // [A * m / kg]
    description = "Amps to weight/length ratio for the cable.";
  }
  
  ComplexNumber dV, Vo;
  ComplexNumber Z;
  ComplexNumber I;
  
  int n_phases {
    value = 3; IOstatus = "input"; units = "none";
    description = "Number of electrical phase currents carried by this cable (cable assembly).";
  }
  
  int n_parallel {
    value = 6; IOstatus = "input";units = "none";
    description = "Constraint: max allowed number of parallel wires to make the cable.";
  }
  
  real drop_tolerance {
    value = 5.; IOstatus = "input"; units = "none";
    description = "Constraint: max voltage allowed to drop across cable, as a percentage of input voltage.";
  }
  
  string cable_size {
    value = ""; IOstatus = "output"; units = "none";
    description = "Cable size in AWG or kcmil.";
  }
  
  real paralleled {
    value = 1.; IOstatus = "output"; units = "none";
    description = "Number of wires paralleled to make this cable.";
  }
  
  real ampacity {
    value = 1000.; IOstatus = "output"; units = "A"; // [amps]
    description = "Current carrying capacity for this cable solution.";
  }
  
  int params_test {
    value = 0; IOstatus = "output"; units = "none";
    description = "Error flag for the sizing tool.";
  }
  
  //------------------------------------
  // ****** OPTION VARIABLE SETUP ******
  //------------------------------------

  Option switchDes {
    allowedValues = { "DESIGN", "OFFDESIGN" }
    description = "Determines if the element is in design or off-design mode.";
    rewritableValues = FALSE;
    trigger = TRUE;
  }

  Option switchThermPort {
    allowedValues = { "TRUE", "FALSE" }
    description = "Determines if component needs thermal port.";
    rewritableValues = FALSE; // enables converter optimization
    trigger = TRUE;
  }
  
  
  //----------------------------------------------------------
  // ****** SETUP PORTS, FLOW STATIONS, SOCKETS, TABLES ******
  //----------------------------------------------------------  

  /* ELECTRICAL PORTS */

  ElectricInputPort EP_I {
    description = "Electric input port.";
  }

  ElectricOutputPort EP_O {
    description = "Electric output port.";
  }

  /* SOCKETS */

  Socket S_eThermMass {
    allowedValues = { "Q_heat", "Mass" }
    description = "Thermal mass socket.";
    socketType = "EThermalMass";
  }

  //-------------------------------------------
  // ****** VARIABLE CHANGED METHODOLOGY ******
  //-------------------------------------------

  void variableChanged(string name, any oldVal) {
    if (name == "switchThermPort") {
      if (switchThermPort == "TRUE") {
        create("", "ThermalInputPort", "Q_I");
      }
    } else if (name =="switchDes") {
      if (switchDes == "DESIGN") {
        R.IOstatus = "input";
        L.IOstatus = "input";
        length.IOstatus = "input";
        SpecificAmperage.IOstatus = "input";
      } else if (switchDes == "OFFDESIGN") {
        R.IOstatus = "output";
        L.IOstatus = "output";
        length.IOstatus = "output";
        SpecificAmperage.IOstatus = "output";
      }
    }
  }
  
  //----------------------------------------
  // ****** CABLE SIZING CALCULATIONS ******
  //----------------------------------------
  
  // AC cable metrics function
  void ac_cable_metrics() {

    // cable sizes available
    string cable_sizes[] = { "18 Gauge",  "16 Gauge",  "14 Gauge",   "12 Gauge",
                             "10 Gauge",   "8 Gauge",   "6 Gauge",    "4 Gauge",
                              "2 Gauge",   "1 Gauge",
                                  "1/0",       "2/0",       "3/0",        "4/0",
                            "262 kcmil", "313 kcmil", "373 kcmil",  "444 kcmil",
                            "535 kcmil", "646 kcmil", "777 kcmil", "1111 kcmil" };
      
    // cable mass [kg/m]
    real cable_mass[] = { 0.047621, 0.050598, 0.065479, 0.078873,
                          0.101195, 0.142864, 0.193461, 0.312514,
                          0.467283, 0.584848, 0.721759, 0.886945,
                          1.055107, 1.244104, 1.510485, 1.753056, 
                          2.083428, 2.461421, 2.983766, 3.494206,
                          4.031433, 5.900565 };
               
    // cable reactance [ohms/m]
    real cable_reactances[] = { 0.000213, 0.000203, 0.000187, 0.000174,
                                0.000164, 0.000157, 0.000148, 0.000128,
                                0.000121, 0.000118, 0.000115, 0.000112,
                                0.000112, 0.000108, 0.000112, 0.000108,
                                0.000105, 0.000102, 0.000102, 0.000102,
                                0.000098, 0.000098 };

    // 400 Hz characteristics -> Resistance [ohms/m]
    real cable_ampacities_400[] = {  17,  25,  40,  48, 
                                     62,  77, 103, 136,
                                    179, 204, 236, 268,
                                    299, 337, 377, 406,
                                    443, 485, 522, 557,
                                    614, 698 };

    real cable_resistivities_400[] = { 0.023617, 0.014807, 0.009336, 0.005866,
                                       0.003704, 0.002278, 0.001433, 0.000975,
                                       0.000606, 0.000490, 0.000402, 0.000333,
                                       0.000280, 0.000241, 0.000216, 0.000195,
                                       0.000177, 0.000167, 0.000156, 0.000143,
                                       0.000129, 0.000108 };

    // 680 Hz characteristics
    real cable_ampacities_680[] = {  17,  25,  40, 48,
                                     62,  77, 102, 135,
                                    175, 198, 225, 252,
                                    276, 306, 340, 363,
                                    395, 432, 466, 498,
                                    550,627 };

    real cable_resistivities_680[] = { 0.023617, 0.014808, 0.009338, 0.005867,
                                       0.003708, 0.002283, 0.001444, 0.000991,
                                       0.000631, 0.000522, 0.000438, 0.000376,
                                       0.000328, 0.000292, 0.000266, 0.000244,
                                       0.000223, 0.000211, 0.000196, 0.000179,
                                       0.000161, 0.000134 };

    // 800 Hz characteristics
    real cable_ampacities_800[] = { 17,  25,  40,  48,
                                    62,  77, 102, 135,
                                   173, 195, 221, 246,
                                   268, 296, 328, 350,
                                   381, 418, 451, 482,
                                   532, 606 };

    real cable_resistivities_800[] = { 0.023617, 0.014809, 0.009339, 0.005868,
                                       0.003710, 0.002288, 0.001449, 0.000999,
                                       0.000645, 0.000538, 0.000457, 0.000396,
                                       0.000347, 0.000312, 0.000285, 0.000262,
                                       0.000239, 0.000225, 0.000209, 0.000192,
                                       0.000173, 0.000143 };
    
    // combined characteristics
    real cable_ampacities[][];
    cable_ampacities.append(cable_ampacities_400);
    cable_ampacities.append(cable_ampacities_680);
    cable_ampacities.append(cable_ampacities_800);

    real cable_resistivities[][];
    cable_resistivities.append(cable_resistivities_400);
    cable_resistivities.append(cable_resistivities_680);
    cable_resistivities.append(cable_resistivities_800);

    // perform error check on input data
    if (frequency != 400) {
      if (frequency != 680) {
        if (frequency != 800) {
          cerr << "Incompatible frequency: "
               << frequency << endl
               << "Voltage in: " << EP_I.V.mag << endl
               << "Voltage out: " << EP_O.V.mag << endl;
        }
      }
    }

    // determine number of lines
    real lines;
    if (n_phases == 1) {
      lines = 2;
    } else if (n_phases == 3) {
      lines = 3;
    } else if (n_phases == 9) {
      lines = 9;
    }
    
    // loop through searching for the cable
    params_test = 0; // error check (ie. no cable is found)
    cable_size = "";
    paralleled = 0;
    Mass = inf;
    real temp_current = 0; // temporary test current
    int found;
    int optimal_wire = -1;
    int conductors = 1;
    int search_stop = 0;
    int index = 0;
    int row;

    real percent;
    real resistance; // [ohms/m]
    real reactance; // [ohms/m]
    ComplexNumber impedance; // [ohms]
    ComplexNumber voltage_drop;
    ComplexNumber voltage;
    ComplexNumber current;
    voltage.copy("EP_I.VLN");
    
    if (frequency == 400) {
      row = 0;
    } else if (frequency == 680) {
      row = 1;
    } else if (frequency == 800) {
      row = 2;
    }
    
    while (search_stop == 0) {
      found = 0;
      while (found == 0) {
        if (conductors > n_parallel) {
          search_stop = 1;
          found = 1;
          if (optimal_wire < 0) {
            params_test = 1;
            Mass = 0;
            impedance.setrj(0, 0);
            voltage_drop.setrj(0, 0);
            ampacity = 0;
            cable_size = "";
            cerr << "\n[ERROR]: No cables were found for these settings" 
                 << "\n         Component:              " << getPathName()
                 << "\n         LL voltage:             " << voltage.mag
                 << "\n         Current:                " << current.mag
                 << "\n         Parallelled Conductors: " << n_parallel
                 << "\n         Voltage Drop Tolerance: " << drop_tolerance
                 << endl;
          } 
        } else if (index > cable_sizes.entries()-1) {
          // increase the conductors and start searching again
          conductors++;
          index = 1;
        } else {
          temp_current = cable_ampacities[row][index] * conductors;
          
          // exit loop once compatible cable setup is found
          found = 1;

          resistance = cable_resistivities[row][index]; // [ohms/m]
          reactance = (frequency/60) * cable_reactances[index]; // [ohms/m]
          impedance.setrj(resistance, reactance);
          impedance.setMagPhase(impedance.mag * length / conductors, impedance.phase); // [ohms]
          current.copy("EP_I.I");
          
          percent = 100. * voltage_drop.mag / voltage.mag;
            
          if (temp_current > current.mag) {
            // check to see if cable setup is optimal
            // this is an AC bus, so there may or may not be send & return
            // lines multiplied by the number of paralleled conductors
            /*
              if (percent > drop_tolerance) {
                index = index + 1;
              }
            */
            if ((lines * cable_mass[index] * conductors * length) < Mass) {
              // if optimal weight, then re-initialize weight, index, and
              // number of paralleled conductors
              Mass = lines * cable_mass[index] * conductors * length; // [lbs]
              optimal_wire = index;
              paralleled = conductors;
              index = cable_sizes.entries();  // terminate loop
            } else {
              index = cable_sizes.entries();  // terminate loop
            }
          } else {
            index++; // terminate loop
          }
        }
      }
    }
    
    if (optimal_wire >= 0) {
      // Initialize rest of output
      cable_size = cable_sizes[optimal_wire];
      
      // Single wire resistance and reactance
      resistance = cable_resistivities[row][optimal_wire]; // [ohms/m]
      reactance = (frequency/60) * cable_reactances[optimal_wire]; // [ohms/m]
      
      // Bundled impedance
      impedance.setrj(resistance, reactance);
      impedance.setMagPhase(impedance.mag * length / paralleled, impedance.phase); // [ohms]
      voltage_drop.copy("current");
      voltage_drop.times("impedance");
      Mass = lines * cable_mass[optimal_wire] * length * paralleled;
      ampacity = cable_ampacities[row][optimal_wire] * paralleled;
      
      R = resistance;
      L = reactance / (2*PI*frequency);
    }
  }
  
  // DC cable metrics function
  void dc_cable_metrics() {

    // cable sizes available
    string cable_sizes[] = {  "18 Gauge",  "16 Gauge",  "14 Gauge",   "12 Gauge",
                              "10 Gauge",   "8 Gauge",   "6 Gauge",    "4 Gauge",
                               "2 Gauge",   "1 Gauge",
                                   "1/0",       "2/0",       "3/0",        "4/0",
                             "262 kcmil", "313 kcmil", "373 kcmil",  "444 kcmil",
                             "535 kcmil", "646 kcmil", "777 kcmil", "1111 kcmil" };
      
    // cable resistance [ohms/m]
    real cable_resistivities[] = { 0.023617, 0.014807, 0.009335, 0.005864,
                                   0.003701, 0.002274, 0.001427, 0.000965,
                                   0.000590, 0.000472, 0.000377, 0.000302,
                                   0.000242, 0.000194, 0.000162, 0.000135,
                                   0.000114, 0.000099, 0.000085, 0.000073,
                                   0.000063, 0.000044 };
    // cable Mass [kg/m]
    real cable_mass[] = { 0.047621, 0.050598, 0.065479, 0.078873,
                          0.101195, 0.142864, 0.193461, 0.312514,
                          0.467283, 0.584848, 0.721759, 0.886945,
                          1.055107, 1.244104, 1.510485, 1.753056,
                          2.083428, 2.461421, 2.983766, 3.494206,
                          4.031433, 5.900565 };

    // cable ampacity [A]
    real cable_ampacities[] = {  17,  25,  40,  48,
                                 62,  77, 103, 137,
                                181, 208, 243, 281,
                                321, 376, 436, 487,
                                553, 630, 709, 783,
                                881, 1098 };
    
    // loop through searching for the cable
    params_test = 0; // error check (ie. no cable is found)
    Mass = inf;
    real resistance = 0;
    real voltage_drop = 0;
    cable_size = "";
    int optimal_wire = 0;
    paralleled = 0;
    int conductors = 1;
    real temp_current = 0; // temporary test current
    int search_stop = 0;
    int found;
    int index = 1;
    
    real voltage = EP_I.V.r;
    real current = EP_I.I.r;
    real percent;
          
    while (search_stop == 0) {
      found = 0;
      while (found == 0) {
        if (conductors > n_parallel) {
          search_stop = 1;
          found = 1;
          if (optimal_wire == 0) {
            params_test = 1;
            Mass = 0;
            resistance = 0;
            cable_size = "";
            cerr << "\n[ERROR]: No cables were found for these settings" 
                 << "\n         Component:              " << getPathName()
                 << "\n         LL voltage:             " << voltage
                 << "\n         Current:                " << current * 1000
                 << "\n         Parallelled Conductors: " << n_parallel
                 << "\n         Voltage Drop Tolerance: " << drop_tolerance
                 << endl;
          }
        } else if (index > cable_sizes.entries()) {
          // increase the conductors and start searching again
          conductors++;
          index = 1;
        } else {
          temp_current = cable_ampacities[index-1] * conductors;
          if (temp_current > current) {

            // exit loop once compatible cable setup is found
            found = 1;
            
            // check to see if cable setup is optimal

            // this is a DC bus so there will be send & return lines multiplied
            // by the number of paralleled conductors
            voltage_drop = current * (2 * cable_resistivities[index-1] * length / conductors); // [ohms]
             
            percent = 100 * voltage_drop / voltage;
            
            if (percent > drop_tolerance) {
              index++;
            } else if ((cable_mass[index-1] * conductors * 2) < Mass) {
              // if optimal weight, then re-initialize weight, index, and
              // number of paralleled conductors
              Mass = cable_mass[index-1] * conductors * 2;
              optimal_wire = index;
              paralleled = conductors;
              index = cable_sizes.entries() + 1;  // terminate loop
            }
            else {
              index = cable_sizes.entries() + 1;  // terminate loop
            }
          }
          else {
            index++;  // terminate loop
          }
        }
      }
    }
    // initialize rest of output
    if (optimal_wire != 0) {
      resistance = 2 * cable_resistivities[optimal_wire-1] * length / paralleled; // [ohms]
      cable_size = cable_sizes[optimal_wire-1];
      voltage_drop = current * (2 * cable_resistivities[optimal_wire-1] * length / paralleled); // [ohms]
      Mass = 2 * cable_mass[optimal_wire-1] * length * paralleled;
      ampacity = cable_ampacities[optimal_wire-1] * paralleled;
      
      R = resistance;
    }   
  }

  //-----------------------------------------------
  // ****** PERFORM ENGINEERING CALCULATIONS ******
  //-----------------------------------------------
  
  void calculate() {
    
    // get frequency from input port
    frequency = EP_I.frequency;
    
    // If on-design mode, must select the cable
    if (switchDes == "DESIGN") {
      if ((EP_I.ElectricPowerType == "DC") && (EP_O.ElectricPowerType == "DC")) {
        dc_cable_metrics();
      } else {
        ac_cable_metrics();
      }
    }
    
    // dV = I * Z
    Z.setrj(R,  L * (2*PI*EP_I.frequency));
    dV.copy("EP_I.I");
    dV.times("Z");

    // Vo_ln = Vi_ln - dV_ln
    Vo.copy("EP_I.VLN");
    Vo.minus("dV");

    if ((EP_I.ElectricPowerType == "AC3") && (EP_O.ElectricPowerType == "AC3")) {
      // Convert this VLN we have calculated into a VLL (port expects VLL inputs)
      Vo.ln2ll();
    } // If not AC3, then we do not need to convert, because VLL = VLN

    // set the conditions in the port
    EP_O.frequency = EP_I.frequency;
    EP_O.setIVRMS(EP_I.I.r, EP_I.I.j, Vo.r, Vo.j);
    
    Loss_r = EP_I.S.r - EP_O.S.r;
    Loss_j = EP_I.S.j - EP_O.S.j;
    
    eff = (EP_O.S.r / EP_I.S.r)**sign(EP_I.S.r);	

    real KW_PER_BTU_PER_SEC = 1.05505585;
    Q_heat = sqrt(Loss_r**2 + Loss_j**2);
    Q_heat /= KW_PER_BTU_PER_SEC;

    if (switchThermPort == "TRUE") {
      // run the thermal mass model
      if (!S_eThermMass.isEmpty()) {
        S_eThermMass.execute();
      }
    }

    // TODO: Run this calculation past Dave
    SpecificAmperage = ampacity * length / Mass;
  }
}
#endif