Commit 95f04c4b authored by Blais, Chris's avatar Blais, Chris
Browse files

add minimal layout

parent d8302ae1
Loading
Loading
Loading
Loading

fig/minimal_layout.png

0 → 100644
+55.8 KiB
Loading image diff...
+175 −0
Original line number Diff line number Diff line
%% Cell type:markdown id:0860494d tags:

#

%% Cell type:code id:3db1995c tags:

``` python
import sympy as sy
import numpy as np
import os
import sys
import scipy as sp
import copy
import yaml
from PIL import Image
from constraint import Problem

from itertools import product

from sympy.core.numbers import int_valued

module_path = os.path.join(os.path.dirname(os.path.abspath('')))
if module_path not in sys.path:
    sys.path.insert(0, module_path)
else:
    print("path already in sys.path")

import matplotlib.pyplot as plt
from numerical_labelling.labelling import *
```

%% Output

    path already in sys.path

%% Cell type:code id:157a872c tags:

``` python
from numerical_labelling.constraints import *
```

%% Cell type:code id:efd7d943 tags:

``` python
from math import gcd

def lcm(a, b):
    return abs(a*b) // gcd(a, b)

def get_consts(cfa, cfb):
    lcm_value = lcm(cfa, cfb)
    const_1 = lcm_value / cfa
    const_2 = lcm_value / cfb

    return const_1, const_2
```

%% Cell type:code id:4f5864fc tags:

``` python
def get_varbs_in_eq(equation_list, all_constr):
    """
    equation_list: list of ids in all constr you want to use
    """
    varb_list = []
    for iconstr in equation_list:
        constr = all_constr[iconstr]
        symbs = [symbol.name for symbol in list(constr.free_symbols)]
        for symb in symbs:
            if symb not in varb_list:
                varb_list.append(symb)
    varb_list = sorted(varb_list)
    return varb_list
```

%% Cell type:markdown id:10524321 tags:

## Use minimal example:
<img src="../fig/minimal_layout.png" width="75%">

%% Cell type:code id:bdd8dd2b tags:

``` python
# just adapt ccro operating point
x0_ccro_file = os.path.join("./saved_x0", "x0_ccro_final_ans_rev2.yaml")

with open(x0_ccro_file, "r") as f:
    x0_ccro = yaml.load(f, yaml.FullLoader)

# convert from ccro to minimal
# cannot use pressure in this example, doesn't balance without pump
x0_min = {}
x0_min["m1"] = x0_ccro["m2"]
x0_min["x1"] = x0_ccro["x2"]
# x0_min["p1"] = x0_ccro["p2"]

x0_min["m2"] = x0_ccro["m3"]
x0_min["x2"] = x0_ccro["x3"]
# x0_min["p2"] = x0_ccro["p3"]

x0_min["m3"] = x0_ccro["m4"]
x0_min["x3"] = x0_ccro["x4"]
# x0_min["p3"] = x0_ccro["p4"]

x0_min["m4"] = x0_ccro["m5"]
x0_min["x4"] = x0_ccro["x5"]
# x0_min["p4"] = x0_ccro["p5"]

x0_min["m5"] = x0_ccro["m6"]
x0_min["x5"] = x0_ccro["x6"]
# x0_min["p5"] = x0_ccro["p6"]

x0_min["m6"] = x0_ccro["m8"]
x0_min["x6"] = x0_ccro["x8"]
# x0_min["p6"] = x0_ccro["p8"]
```

%% Cell type:code id:0895969d tags:

``` python
# check
m1 = x0_min["m1"]
m2 = x0_min["m2"]
m3 = x0_min["m3"]
m4 = x0_min["m4"]
m5 = x0_min["m5"]
m6 = x0_min["m6"]


x1 = x0_min["x1"]
x2 = x0_min["x2"]
x3 = x0_min["x3"]
x4 = x0_min["x4"]
x5 = x0_min["x5"]
x6 = x0_min["x6"]

# check that constraints equal 0
print (m1 - m2 + m6)
print(m1*x1 - m2*x2 + m6*x6)

print (m2 - m3 - m4)
print(m2*x2 - m3*x3 - m4*x4)

print (m4 - m5 - m6)
print(m4*x4 - m5*x5 - m6*x6)
```

%% Output

    0.0
    0.0
    0.0
    0.0
    0.0
    0.0

%% Cell type:code id:793bf1e7 tags:

``` python
os.path.exists("../testing")
```

%% Output

    True

%% Cell type:code id:a08e9e01 tags:

``` python
# save
x0_min_file = os.path.join("../testing/operating_points/", "x0_minimal.yaml")

with open(x0_min_file, "w") as f:
    yaml.dump(x0_min, f)
```
+12 −0
Original line number Diff line number Diff line
m1: 50.0
m2: 70.0
m3: 10.0
m4: 60.0
m5: 40.0
m6: 20.0
x1: 43.0
x2: 45.0
x3: 15.0
x4: 50.0
x5: 50.0
x6: 50.0
+328 −0
Original line number Diff line number Diff line
import sympy as sy
import os
import sys

module_path = os.path.join(os.path.dirname(os.path.abspath('')))
if module_path not in sys.path:
    sys.path.insert(0, module_path)
else:
    print("path already in sys.path")

from numerical_labelling.constraints import *

def constr_minimal():
    """
    minimal example: membrane, 6 streams, recycle
    returns
    constr: list of symbolic constraints
    varbs_sys: list of the process variables (flow, concentration, pressure)
    coeff_sys: list of the constant coefficients
    """

    m1, m2, m3, m4, m5, m6 = sy.symbols(
        "m1 m2 m3 m4 m5 m6"
    )
    x1, x2, x3, x4, x5, x6= sy.symbols(
        "x1 x2 x3 x4 x5 x6"
    )
    p1, p2, p3, p4, p5, p6 = sy.symbols(
        "p1 p2 p3 p4 p5 p6"
    )
    varbs_list_sys = [
        m1, m2, m3, m4, m5, m6,
        x1, x2, x3, x4, x5, x6,
        p1, p2, p3, p4, p5, p6,
    ]
    varbs_sys = {varb.name:varb for varb in varbs_list_sys}

    # membrane constants
    a1, a1p, a1pp = sy.symbols(
        'a1 a1p a1pp'
    )
    b1, b1p = sy.symbols(
        'b1 b1p'
    )
    c1, c1p = sy.symbols(
        'c1 c1p'
    )

    coeff_sys_list = [
        a1, a1p, a1pp,
        b1, b1p,
        c1, c1p,
    ]

    coeff_sys = {coeff.name:coeff for coeff in coeff_sys_list}

    # setup pilot constraints
    constr_sys = {}
    sltns_sys = {}

    constr_sys, sltns_sys = constr_mix(
        min1=m1, min2=m6, mout=m2,
        xin1=x1, xin2=x6, xout=x2,
        pin1=p1, pin2=p6, pout=p2,
        constr=constr_sys, sltns=sltns_sys,
    )


    # membrane 1
    constr_sys, sltns_sys = constr_mex(
        mfeed=m2, mperm=m3, mret=m4,
        xfeed=x2, xperm=x3, xret=x4,
        pfeed=p2, pperm=p3, pret=p4,
        a=a1, ap=a1p, app=a1pp,
        b=b1, bp=b1p,
        c=c1, cp=c1p,
        constr=constr_sys, sltns=sltns_sys,
    )

    # split point
    constr_sys, sltns_sys = constr_split(
        min=m4, mout1=m5, mout2=m6,
        xin=x4, xout1=x5, xout2=x6,
        pin=p4, pout1=p5, pout2=p6,
        constr=constr_sys, sltns=sltns_sys
    )

    return constr_sys, varbs_sys, coeff_sys

def constr_ccro():
    """
    returns
    constr: list of symbolic constraints
    varbs_sys: list of the process variables (flow, concentration, pressure)
    coeff_sys: list of the constant coefficients
    """

    m1, m2, m3, m4, m5, m6, m7, m8 = sy.symbols(
        "m1 m2 m3 m4 m5 m6 m7 m8"
    )
    x1, x2, x3, x4, x5, x6, x7, x8 = sy.symbols(
        "x1 x2 x3 x4 x5 x6 x7 x8"
    )
    p1, p2, p3, p4, p5, p6, p7, p8 = sy.symbols(
        "p1 p2 p3 p4 p5 p6 p7 p8"
    )
    varbs_list_sys = [
        m1, m2, m3, m4, m5, m6, m7, m8,
        x1, x2, x3, x4, x5, x6, x7, x8,
        p1, p2, p3, p4, p5, p6, p7, p8,
    ]
    varbs_sys = {varb.name:varb for varb in varbs_list_sys}

    # membrane constants
    a1, a1p, a1pp = sy.symbols(
        'a1 a1p a1pp'
    )
    b1, b1p = sy.symbols(
        'b1 b1p'
    )
    c1, c1p = sy.symbols(
        'c1 c1p'
    )

    # pump constants
    e1, h1, g1, j1, e2, h2, g2, j2 = sy.symbols("e1 h1 g1 j1 e2 h2 g2 j2")

    coeff_sys_list = [
        a1, a1p, a1pp,
        b1, b1p,
        c1, c1p,
        e1, h1, g1, j1, e2, h2, g2, j2,
    ]

    coeff_sys = {coeff.name:coeff for coeff in coeff_sys_list}

    # setup pilot constraints
    constr_sys = {}
    sltns_sys = {}

    constr_sys, sltns_sys = constr_pump(
        min = m1, mout=m2, 
        xin=x1, xout = x2,
        pin=p1, pout=p2,
        E=e1, H = h1, G = g1, J= j1,
        constr = constr_sys, sltns=sltns_sys,
    )

    constr_sys, sltns_sys = constr_mix(
        min1=m2, min2=m8, mout=m3,
        xin1=x2, xin2=x8, xout=x3,
        pin1=p2, pin2=p8, pout=p3,
        constr=constr_sys, sltns=sltns_sys,
    )


    # membrane 1
    constr_sys, sltns_sys = constr_mex(
        mfeed=m3, mperm=m4, mret=m5,
        xfeed=x3, xperm=x4, xret=x5,
        pfeed=p3, pperm=p4, pret=p5,
        a=a1, ap=a1p, app=a1pp,
        b=b1, bp=b1p,
        c=c1, cp=c1p,
        constr=constr_sys, sltns=sltns_sys,
    )

    # split point
    constr_sys, sltns_sys = constr_split(
        min=m5, mout1=m6, mout2=m7,
        xin=x5, xout1=x6, xout2=x7,
        pin=p5, pout1=p6, pout2=p7,
        constr=constr_sys, sltns=sltns_sys
    )

    constr_sys, sltns_sys = constr_pump(
        min=m7, mout=m8, 
        xin=x7, xout = x8,
        pin=p7, pout=p8,
        E=e2, H = h2, G = g2, J= j2,
        constr = constr_sys, sltns=sltns_sys,
    )

    return constr_sys, varbs_sys, coeff_sys


def constr_pilot():
    """
    returns
    constr: list of symbolic constraints
    varbs_sys: list of the process variables (flow, concentration, pressure)
    coeff_sys: list of the constant coefficients
    """
    m1, m2, m3, m4, m5, m6, m7, m8, m9, m10, m11, m12, m13, m14, m15, m16, m17 = sy.symbols(
        "m1 m2 m3 m4 m5 m6 m7 m8 m9 m10 m11 m12 m13 m14 m15 m16 m17"
    )
    x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x13, x14, x15, x16, x17 = sy.symbols(
        "x1 x2 x3 x4 x5 x6 x7 x8 x9 x10 x11 x12 x13 x14 x15 x16 x17"
    )
    p1, p2, p3, p4, p5, p6, p7, p8, p9, p10, p11, p12, p13, p14, p15, p16, p17 = sy.symbols(
        "p1 p2 p3 p4 p5 p6 p7 p8 p9 p10 p11 p12 p13 p14 p15 p16 p17"
    )
    varbs_list_sys = [
        m1, m2, m3, m4, m5, m6, m7, m8, m9, m10, m11, m12, m13, m14, m15, m16, m17,
        x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x13, x14, x15, x16, x17,
        p1, p2, p3, p4, p5, p6, p7, p8, p9, p10, p11, p12, p13, p14, p15, p16, p17,
    ]
    varbs_sys = {varb.name:varb for varb in varbs_list_sys}

    # membrane constants
    a1, a1p, a1pp, a2, a2p, a2pp, a3, a3p, a3pp, a4, a4p, a4pp = sy.symbols(
        'a1 a1p a1pp a2 a2p a2pp a3 a3p a3pp a4 a4p a4pp'
    )
    b1, b1p, b2, b2p, b3, b3p, b4, b4p = sy.symbols(
        'b1 b1p, b2 b2p b3 b3p b4 b4p'
    )
    c1, c1p, c2, c2p, c3, c3p, c4, c4p = sy.symbols(
        'c1 c1p c2 c2p c3 c3p c4 c4p'
    )

    # pump constants
    e1, h1, g1, j1, k1, l1 = sy.symbols("e1 h1 g1 j1 k1 l1")

    coeff_sys_list = [
        a1, a1p, a1pp, a2, a2p, a2pp, a3, a3p, a3pp, a4, a4p, a4pp,
        b1, b1p, b2, b2p, b3, b3p, b4, b4p,
        c1, c1p, c2, c2p, c3, c3p, c4, c4p,
        e1, h1, g1, j1, k1, l1
    ]

    coeff_sys = {coeff.name:coeff for coeff in coeff_sys_list}

    # setup pilot constraints
    constr_sys = {}
    sltns_sys = {}
    constr_sys, sltns_sys = constr_mix(
        min1=m1, min2=m17, mout=m2,
        xin1=x1, xin2=x17, xout=x2,
        pin1=p1, pin2=p17, pout=p2,
        constr=constr_sys, sltns=sltns_sys,
    )

    constr_sys, sltns_sys = constr_pump(
        min = m2, mout=m3, 
        xin=x2, xout = x3,
        pin=p2, pout=p3,
        E=e1, H = h1, G = g1, J= j1,
        constr = constr_sys, sltns=sltns_sys,
    )
    # membrane 1
    constr_sys, sltns_sys = constr_mex(
        mfeed=m3, mperm=m4, mret=m5,
        xfeed=x3, xperm=x4, xret=x5,
        pfeed=p3, pperm=p4, pret=p5,
        a=a1, ap=a1p, app=a1pp,
        b=b1, bp=b1p,
        c=c1, cp=c1p,
        constr=constr_sys, sltns=sltns_sys,
    )
    # membrane 2
    constr_sys, sltns_sys = constr_mex(
        mfeed=m5, mperm=m6, mret=m8,
        xfeed=x5, xperm=x6, xret=x8,
        pfeed=p5, pperm=p6, pret=p8,
        a=a2, ap=a2p, app=a2pp,
        b=b2, bp=b2p,
        c=c2, cp=c2p,
        constr=constr_sys, sltns=sltns_sys,
    )
    # membrane 3
    constr_sys, sltns_sys = constr_mex(
        mfeed=m8, mperm=m9, mret=m11,
        xfeed=x8, xperm=x9, xret=x11,
        pfeed=p8, pperm=p9, pret=p11,
        a=a3, ap=a3p, app=a3pp,
        b=b3, bp=b3p,
        c=c3, cp=c3p,
        constr=constr_sys, sltns=sltns_sys,
    )
    # membrane 4
    constr_sys, sltns_sys = constr_mex(
        mfeed=m11, mperm=m12, mret=m14,
        xfeed=x11, xperm=x12, xret=x14,
        pfeed=p11, pperm=p12, pret=p14,
        a=a4, ap=a4p, app=a4pp,
        b=b4, bp=b4p,
        c=c4, cp=c4p,
        constr=constr_sys, sltns=sltns_sys,
    )

    # mixing points
    constr_sys, sltns_sys = constr_mix(
        min1=m4, min2=m6, mout=m7,
        xin1=x4, xin2=x6, xout=x7,
        pin1=p4, pin2=p6, pout=p7,
        constr=constr_sys, sltns=sltns_sys
    )
    constr_sys, sltns_sys = constr_mix(
        min1=m7, min2=m9, mout=m10,
        xin1=x7, xin2=x9, xout=x10,
        pin1=p7, pin2=p9, pout=p10,
        constr=constr_sys, sltns=sltns_sys
    )
    constr_sys, sltns_sys = constr_mix(
        min1=m10, min2=m12, mout=m13,
        xin1=x10, xin2=x12, xout=x13,
        pin1=p10, pin2=p12, pout=p13,
        constr=constr_sys, sltns=sltns_sys
    )

    # split point
    constr_sys, sltns_sys = constr_split(
        min=m14, mout1=m15, mout2=m16,
        xin=x14, xout1=x15, xout2=x16,
        pin=p14, pout1=p15, pout2=p16,
        constr=constr_sys, sltns=sltns_sys
    )

    constr_sys, sltns_sys = constr_valve(
        min=m16, mout=m17,
        xin=x16, xout=x17,
        pin=p16, pout=p17,
        K=k1,L=l1,
        constr=constr_sys, sltns=sltns_sys,
    )
    
    return constr_sys, varbs_sys, coeff_sys
+24 −0
Original line number Diff line number Diff line
@@ -22,6 +22,30 @@ def build_pgraph(sys_name="ccro", overwrite=False, n_components=3, arc_splitter=
        else:
            raise Exception(f"ERROR: file {pgraph_path} already exists and overwrite not specified. Aborting.")

    if sys_name == "minimal":
        # almost the exact same layout as CCRO, but combining streams
        incidence = np.array(
            [   # e  Mx  Me1 Sp1  P2
                [-1, +1,  0,  0],  # 1
                [ 0, -1, +1,  0],  # 2
                [+1,  0, -1,  0],  # 3
                [ 0,  0, -1, +1],  # 4
                [+1,  0,  0, -1],  # 5
                [ 0, +1,  0, -1],  # 6
            ]
        )

        coordinates = np.array([
            [2, 1], #env
            [0.75, 0], #mx
            [1.5, 0], #me1
            [2, -1],] # sp1
        )
        if arc_splitter is None:
            arc_splitter =[[4 - 1, 5 - 1, 6 - 1]]
        if arc_equal_composition is None:
            arc_equal_composition = [[4 - 1, 5 - 1, 6 - 1]]

    if sys_name == "ccro":
        incidence = np.array(
            [   # e   P1  Mx Me1 Sp1  P2
Loading