Simple circuit

Preliminary note:

This case is taken from OpenMDAO. OpenMDAO is an open-source computing platform for systems analysis and multidisciplinary optimization developed by the NASA. Its philosophy shares some of the goals with CoSApp. So this example is also an opportunity to compare both software.

OpenMDAO is licensed under Apache License.

Case description

This is a simple circuit having a source providing the circuit with a constant intensity. The resolution is driving the potential V at node1 and node2 so that the current flowing in equals the one flowing out.

simple-circuit

Creating elementary bricks

First of all, the 5 elements need to be modeled:

  • Two resistors

  • One diode

  • Two nodes

  • The source

  • The ground

[1]:
from cosapp.base import System, Port
from cosapp.drivers import NonLinearSolver, NonLinearMethods
import numpy as np

class Voltage(Port):
    def setup(self):
        self.add_variable('V')


class Intensity(Port):
    def setup(self):
        self.add_variable('I')


class Resistor(System):

    def setup(self, R=1.):
        self.add_input(Voltage, 'V_in')
        self.add_input(Voltage, 'V_out')
        self.add_output(Intensity, 'I')

        self.add_inward('R', R, desc='Resistance in Ohms')
        self.add_outward('deltaV')  # Not mandatory; could be local to compute method

    def compute(self):
        self.deltaV = self.V_in.V - self.V_out.V
        self.I.I = self.deltaV / self.R


class Diode(System):
    """Diode model

    The current intensity flowing through the diode is calculated based on

    $ I = I_s \\exp \\left( \\dfrac{V_{in} - V_{out}}{V_t} - 1 \\right) $
    """

    def setup(self):
        self.add_input(Voltage, 'V_in')
        self.add_input(Voltage, 'V_out')
        self.add_output(Intensity, 'I')

        self.add_inward('Is', 1e-15, desc='Saturation current in Amps')
        self.add_inward('Vt', .025875, desc='Thermal voltage in Volts')

        self.add_outward('deltaV')  # Not mandatory; could be local to `compute` method

    def compute(self):
        self.deltaV = self.V_in.V - self.V_out.V
        self.I.I = self.Is * np.exp(self.deltaV / self.Vt - 1.)


class Node(System):

    def setup(self, n_in=1, n_out=1):
        self.add_property('n_in', int(n_in))
        self.add_property('n_out', int(n_out))

        if min(self.n_in, self.n_out) < 1:
            raise ValueError("Node needs at least one incoming and one outgoing current")

        self.add_property('incoming', tuple(
            self.add_input(Intensity, f"I_in{i}")
            for i in range(self.n_in)
        ))
        self.add_property('outgoing', tuple(
            self.add_input(Intensity, f"I_out{i}")
            for i in range(self.n_out)
        ))

        self.add_inward('V')
        self.add_outward('sum_I_in', 0., desc='Sum of all incoming currents')
        self.add_outward('sum_I_out', 0., desc='Sum of all outgoing currents')

        # Off-design current balance problem
        self.add_unknown('V')
        self.add_equation('sum_I_in == sum_I_out', name='current balance')

    def compute(self):
        self.sum_I_in = sum(current.I for current in self.incoming)
        self.sum_I_out = sum(current.I for current in self.outgoing)


class Source(System):

    def setup(self, I=0.1):
        self.add_inward('I', I)
        self.add_output(Intensity, 'I_out', {'I': I})

    def compute(self):
        self.I_out.I = self.I


class Ground(System):

    def setup(self, V=0.):
        self.add_inward('V', V)
        self.add_output(Voltage, 'V_out', {'V': V})

    def compute(self):
        self.V_out.V = self.V

Building the circuit

The circuit is built by reusing Node and Resistor elements customized thanks to keywords arguments that will be processed in the setup method of those classes.

[2]:
class Circuit(System):

    def setup(self):
        n1 = self.add_child(Node('n1', n_in=1, n_out=2), pulling={'I_in0': 'I_in'})
        n2 = self.add_child(Node('n2'))

        R1 = self.add_child(Resistor('R1', R=100.), pulling={'V_out': 'Vg'})
        R2 = self.add_child(Resistor('R2', R=10000.))
        D1 = self.add_child(Diode('D1'), pulling={'V_out': 'Vg'})

        self.connect(R1.V_in, n1.inwards, 'V')
        self.connect(R2.V_in, n1.inwards, 'V')
        self.connect(R1.I, n1.I_out0)
        self.connect(R2.I, n1.I_out1)

        self.connect(R2.V_out, n2.inwards, 'V')
        self.connect(D1.V_in, n2.inwards, 'V')
        self.connect(R2.I, n2.I_in0)
        self.connect(D1.I, n2.I_out0)

p = System('model')
p.add_child(Source('source', I=0.1))
p.add_child(Ground('ground', V=0.))
p.add_child(Circuit('circuit'))

p.connect(p.source.I_out, p.circuit.I_in)
p.connect(p.ground.V_out, p.circuit.Vg)

Solving the problem

[3]:
p.add_driver(NonLinearSolver('solver', method=NonLinearMethods.POWELL))

p.run_drivers()

print(p['circuit.n1.V'])
print(p['circuit.n2.V'])
print(p['circuit.R1.I'])
print(p['circuit.R2.I'])
print(p['circuit.D1.I'])

# sanity check: current should sum to 0.1 A
print('sanity check: should sum to 0.1 A: ', p['circuit.R1.I.I'] + p['circuit.D1.I.I'])
9.908302820430945
0.7385848635255148
Intensity: {'I': 0.09908302820430945}
Intensity: {'I': 0.000916971795690543}
Intensity: {'I': 0.0009169717958134371}
sanity check: should sum to 0.1 A:  0.10000000000012288