ABSTRACT (original Article): Mathematica provides unique symbolic and machine learning tools for developing and solving mathematical models that make it attractive for engineering applications. Chemical engineers frequently have a need to update process design software with custom unit operation models. This paper shows two new methods to connect Mathematica to chemical process simulators, including a local PC-based connection and a cloud-based connection. Both methods are discussed in detail using an illustrative example of a membrane separation process. We also demonstrate how the connection between Mathematica and Excel facilitates spreadsheet calculations while enabling the connection with the process simulator software. Furthermore, all methods shown here are model-independent and can be implemented with any equation-based model or with machine-learning models when adequate training data is available. This document shows the Mathematica code development and methods for deployment to Excel. The original article goes further and shows how to connect the code to the chemical process simulator. CITATION (original article): Cowart, S.V., Yuk, S.F., James, C.M., Nagelli, E.A. and Biaglow, A.I. (2024), Augmented and Cloud Computing With Chemical Process Simulators. Softw: Pract Exper. https://doi.org/10.1002/spe.3397

I. Introduction

One way to separate the components in a mixture or solution is through the use of industrial membrane technology. The membrane separation process contains a semipermeable film or barrier that allows selective dissolution and diffusion of the mixture components. The degree of separation is limited by the surface area of the membrane. Engineering calculations are used to determine the surface area of the membrane required for a specified degree of separation. Additionally, though not discussed here, there are power and equipment costs associated with pressurizing the mixture to a level adequate to induce flow through the process. Various strategies are used to model the separation, and the simplest of these is the single-stage perfectly-mixed model (SPM), which is the starting point for the methods described in this paper.
A detailed description of the SPM can be found in the references [1,2]. The model is derived in English units, with units explained in Reference 2. A diagram of the SPM process is shown below. As shown in the figure, two feed streams are fed to the membrane unit, with feed molar flow rates given by
n
F10
and
n
F20
, and feed mole fractions by
z
1
and
z
2
. In this document, feed stream 1 is air containing nitrogen and oxygen, with other species normally found in air, such as water,
CO
2
and trace gases, assumed zero. The molar flow rate of stream 1 is
n
F10
, and
z
1
is a list containing the mole fractions of nitrogen and oxygen in the air, assumed to be 0.79 and 0.21, respectively. Feed stream 2 is a purge stream and the flow rate and mole fractions of stream 2 are
n
F20
and
z
2
, respectively. These values can, in principle, be anything; however, we define
z
2
as a list containing the two components
(
N
2
and
O
2
)
to simplify the cascade model. The membrane itself is shown by the cross-hatching and in this case it is a 0.2 μm thick low-density polyethylene film. An inert porous support is shown by the dotted zone underneath the membrane. The support is necessary to enable the membrane to withstand the moderately high pressures that are applied to the feed, but is considered to be inert with respect to the separation. The gas that goes through the membrane is called the “permeate,” and the gas that does not get through is the “retentate.”
Two mole balance and two molar flux equations describe the transport of nitrogen and oxygen in the membrane (Section II, eq12-eq15). These four equations form a system of four solvable algebraic equations. The approach used here demonstrates a unique feature of Mathematica, namely, the ability to embed an unsolved system of equations into a function and then allow the symbolic or numeric engine to solve the system of equations each time the function is called. In other words, the solution to the system of equations does not need to be found explicitly, and the function can be called over and over again, with the solver embedded into the function code. Furthermore, the function can be called directly from Excel or deployed to the cloud for execution in Excel.
The usefulness of the SPM is that it forms the basis for more advanced configurations such as co-current, cross-current, and countercurrent cascades. The cascade models are created from the SPM model. An example of a two-stage co-current cascade is shown below. The co-current flow pattern refers to the movement of the input and output streams in the same direction.
In the cascade, the individual stages are connected with streams that contain flow rate and composition information. The definition of the single stage as a function allows the function to receive the intermediate stream information, return the output values, and then get re-used in subsequent stages in the cascade, with stream information passed into and out of the single stage function as the solution progresses through the cascade. This document shows the SPM and co-current cascade models.
The rest of this document is organized as follows: Section II shows the code for the SPM and for embedding the model into a function. Section III develops the code to use the SPM function in a co-current multi-stage cascade. As the complexity of the function-based models increases, they can take significant time to run. Data-based models can be much faster, although requiring training data. Accordingly, Sections II-III show the development of the co-current membrane model because this is used as a simulation to generate training data for the machine learning model, and Sections IV-VII show the development of a data-based machine-learning model for the co-current cascade. In Section IV, an input dataset is created using gaussian distributions to describe each of the input variables. In Section V, the input dataset is used to generate an output dataset using the cascade function from Section III as a simulation to generate the output data. Section V also defines the associations between the input and output datasets. In Section VI, a simple neural net is used to model the process and the results are compared to the original function model. In Section VII, the neural net model is then exported for use in other programs. Finally, Section VIII shows the deployment of the various models to Excel for use in CHEMCAD calculations. Data mapping into CHEMCAD is fully explained in the original paper (Reference 2).

II. Program Code for Single Stage

This section contains the code for the single stage function. The inputs to the function include flow rates and compositions from the two input streams that feed the stage. The inputs also include the membrane area and pressures that describe the internal conditions of the stage. The function outputs include flow rates and compositions of the streams leaving the stage.
The function is built with the Mathematica function Module, which allows passage of functional inputs as parameters and definition of local variables, facilitating use of the function in the cascade calculations in the next section. If we were to stop with the single stage, the Module structure would not be necessary [1,2].
The calculation begins with definition of the permeabilities and permeances of oxygen and nitrogen

P
M
A
and
P
M
B
and
P
M
A
and
P
M
B

. These constants are based on Henry’s Law solubility constants and diffusion coefficients, which are available in the literature [3]. Permeance is the permeability divided by the membrane thickness, in this case 0.2 μm or
6.6×
-7
10
ft
. Details of this calculation are shown in Reference 2.
The two molar flux and two mole balance equations for nitrogen and oxygen are defined next. These four equations form a system of solvable algebraic equations. A numeric solution using NSolveValues was used to obtain the solution, but SolveValues could have also been used. While NSolveValues produces an approximate solution, the accuracy is more than adequate and the numerical method saves a little time in the solution compared to using SolveValues.
(*FunctionDefinitionofSingleStageWell-MixedModel*)​​(*A-Oxygen;B-Nitrogen*)​​f2[A2_,nF102_,nF202_,z12_,z22_,PR2_,PP2_]:=​​Module{A=A2,nF10=nF102,nF20=nF202,z1=z12,z2=z22,PR=PR2,PP=PP2,ans},​​​​
P
M
A
=(.472*
-6
10
)*(.46*
-6
10
)*
Times[
]
(*PermeabilityofA*);​​
P
M
B
=(.228*
-6
10
)*(.32*
-6
10
)*
Times[
]
(*PermeabilityofB*);​​
P
M
A
=
P
M
A
(0.66*
-6
10
)(*PermeanceofA*);​​
P
M
B
=
P
M
B
(0.66*
-6
10
)(*PermeanceofB*);​​​​eq12=y*(nF10-nR1)==A*
P
M
A
*(x*PR-y*PP)(*molarfluxofA*);​​eq13=(1-y)*(nF10-nR1)==A*
P
M
B
*((1-x)*PR-(1-y)*PP)(*molarfluxofB*);​​eq14=nF10+nF20==nR1+nP1(*Totalmolebalance*);​​eq15=z1*nF10+z2*nF20==x*nR1+y*nP1(*Speciesmolebalanceon
O
2
*);​​eq16=0<=y<=1(*constraintonmolefractionof
O
2
inpermeate*);​​eqns={eq12,eq13,eq14,eq15,eq16};​​vars={x,y,nP1,nR1};​​NSolveValues[eqns,vars,Reals]
Functions to be used in the co-current staging are derived from f2 as shown below:
In[]:=
x[A2_,nF102_,nF202_,z12_,z22_,PR2_,PP2_]:=f2
[[1]][[1]];​​y[A2_,nF102_,nF202_,z12_,z22_,PR2_,PP2_]:=f2
[[1]][[2]];​​nP1[A2_,nF102_,nF202_,z12_,z22_,PR2_,PP2_]:=f2
[[1]][[3]];​​nR1[A2_,nF102_,nF202_,z12_,z22_,PR2_,PP2_]:=f2
[[1]][[4]];
The definition of the single stage as a function allows the function to be re-used in cascade calculations. Stream information is passed into and out of the single stage function as the solution progresses through the cascade. Accordingly the function developed in this section is used in the staging model of the next section.

III. Program Code for Co-current Staging

Because of the concentration gradients in an actual membrane unit, the SPM model is not truly representative of the actual system unless the membrane area is differentially small. However, as explained in the Introduction, two or more SPM models can be networked into a cascade. The cascade configuration is useful for modeling larger membrane units where the SPM model does not adequately describe the system.
In such a case, the total membrane area can divided by the number of stages in the cascade to obtain the membrane area per stage, with each stage modeled with the SPM. The following diagram illustrates the situation that arises when a single membrane unit is conceptually divided into ten equal sections, and then each section is treated as a single stage with the SPM model in a cascade:
Clearly, the membrane area of each section in such a system, ΔA, will get smaller if the number of sections is increased. Presumably, if more stages are defined and ΔA is further reduced, any concentration gradients would also be reduced further. In the limit as the number of stages approaches infinity, ΔA approaches zero, the concentration gradients will also approach zero, and the resulting differentially small sections could then be exactly described with the perfect mixing model (SPM).
In the development that follows, there are 100 stages, so the area per stage is the total area divided by 100. The number of stages could be increased, but this will increase the calculation time. Likewise, to decrease calculation time, the number of stages can be reduced to a minimum of two, which is useful for testing or proof-of-concept experiments. Or, the differential mole balance and molar flux equations could be integrated. We do not do that here because we have a working SPM model for the discrete stage, the SPM model is simple to encode, and we seek to demonstrate the usefulness of that model. The function f2 from the previous section, which incorporates the SPM model, is used to calculate the flow rates and compositions of the streams leaving each stage, and these outputs then form the inputs for the next stage.
The entire cascade is modeled as a single function f4, built with Module, which takes as inputs the molar flow rate of the air feed, the permeate and retentate pressures, and the total membrane area. These four variables are passed as parameters into the Module. All other streams internal to the cascade are internal to the Module. Also, as with the single stage, the purge feed can be any desired value but is set to zero inside the Module with local variable nF20int.
(*FunctionDefinitionforCocurrentStagingoftheWell-MixedModel*)​​f4[nF_,PR_,PP_,A_]:=​​Module[{nF10int=nF,PRint=PR,PPint=PP,area=A,nF20int,z2int,z1int,size,Alist,xlist,ylist,nP1list,nR1list,xFint,i},​​size=100;(*numberofstages*)​​​​(*Initializelistsforstorageofintermediatestagedata*)​​Alist=NestList[#+1&,Range[1],size-1];(*listofstageareas*)​​xlist=NestList[#+1&,Range[1],size-1];(*retenatemolefractionlist*)​​ylist=NestList[#+1&,Range[1],size-1];(*permeatemolefractionlist*)​​nP1list=NestList[#+1&,Range[1],size-1];(*permeatemolarflowlist*)​​nR1list=NestList[#+1&,Range[1],size-1];(*retenatemolarlowlist*)​​​​(*TaketotalareaanddividebysizeandplacetheresultintoeachpositioninAlist*)​​For[i=1,i<size+1,i++,Alist[[i]]=area/size];​​​​z1int=0.21;(*initalenteringfeedzforSPMfunctionswheni=1*)​​z2int=0.;(*initialenteringpurgezforSPMfunctionswheni=1*)​​nF20int=0.;(*initialenteringpurgeflowforSPMfunctionswheni=1*)​​​​(*IntheForloop,thefunctionsx,y,nR1,andnP1comefromthesinglestagemodelf2*)​​For[i=1,i<size+1,i++,​​xlist[[i]]=Quiet[x[Alist[[i]],nF10int,nF20int,z1int,z2int,PRint,PPint]];​​ylist[[i]]=Quiet[y[Alist[[i]],nF10int,nF20int,z1int,z2int,PRint,PPint]];​​nR1list[[i]]=Quiet[nR1[Alist[[i]],nF10int,nF20int,z1int,z2int,PRint,PPint]];​​nP1list[[i]]=Quiet[nP1[Alist[[i]],nF10int,nF20int,z1int,z2int,PRint,PPint]];​​nF10int=nR1list[[i]];(*re-assignedforpassingtonextstage*)​​nF20int=nP1list[[i]];(*re-assignedforpassingtonextstage*)​​z1int=xlist[[i]];(*re-assignedforpassingtonextstage*)​​z2int=ylist[[i]]];(*re-assignedforpassingtonextstage*)​​​​(*functionoutputs*)​​{(1-(nR1list[[size]]/nF)),xlist[[size]],ylist[[size]],nP1list[[size]],nR1list[[size]]}]
The following function call verifies that everything is working, and can be checked against Reference 3:
In[]:=
f4[3343,150,15,22011.]

IV. Generation of Normally Distributed Input Variables

This section produces input data which is then fed to the functions from the previous section to generate training data for the neural net model. The four input variables include the feed flow rate (nF), the retentate pressure (PR) , the permeate pressure (PP), and the membrane area (A), and these are assigned distributions d1 to d4, respectively. The data for each input variable are randomly generated from a normal distribution having a specified mean and standard deviation. The distributions in the data are “narrow,” having a relative standard deviations of about 1% of the mean. The standard deviation can be increased as desired, but increasing the standard deviation will increase the error in the final neural net model, since it is using a single linear layer.

V. Generation of Output Data and Associations for Training

Once the individual input data lists are created, associated output data lists must be generated. Then, in the structured machine learning paradigm, the associations needed for training the neural net need to be defined. This is done in four steps in the code below. In step 1 the lists are initialized. In step 2, input dataset ds12 is created in a For loop by generating a list of d1 to d4. This is done for each set of four numbers in d1 through d4. So ds12 is a list of lists, where ds12[i]={d1[i], d2[i], d3[i], d4[i]}. Next, in the same For loop, the ds12 values are fed to the model function f4 to generate a list of lists of output data, ds22. Third, the input and output datasets are scaled using appropriate factors to create lists ds13 and ds23. Fourth and finally, the individual lists in ds13 are associated with the corresponding outputs in ds23 and the resulting data are placed into dataset ds30.
The dataset ds30 which is subsequently used for training is viewed below for verification purposes.
This particular calculation only uses 100 of the 10,000 available inputs from the original distributions. This number can be increased up to 10,000, but this takes more time. The 100-row training set took about 10 minutes, while the 10,000-row training set took 5.4 hours, although these times are highly dependent on the computer hardware.

VI. Building, Training, and Results of Neural Net Model

This section builds and trains the neural net model using the training data in ds30. The neural net in this example consists of a single linear layer. The neural net model is defined first using NetChain, followed by training with NetTrain.
The results of the neural net model are then visually compared to the outputs generated by f4. Five different randomly selected input rows from the ds13 input dataset are displayed first. The results from the corresponding output rows from the ds23 output dataset (obtained from function f4) are displayed next. Then the neural net function net is applied to the same five inputs, and the results displayed for point-by-point inspection.
Point-by-point comparison of the last two matrices shows that the results are reasonably close. The next calculations shows the sum-of-squares error and the average sum-of-squares error for the limited results. Accuracy is also calculated using the SSQ error on the same five outputs and then again on the complete dataset, with SSQ being divided by the length of the dataset. The NetTrain command was incorporated in a recursive fashion so that the training can be run multiple times. To obtain these results, the training was conducted three times in succession. SSQ stops decreasing after 3 training runs, and then starts to increase with more training.
Finally, the data are compared using histograms of the calculated output and training output data. The amber-colored bars represent the outputs from f4, the blue-colored bars are the outputs from the neural net, and the gray color represents the overlap.

VII. Saving the Trained Neural Net

Creating the data set and training the neural net can take significant time, especially when large data sets are used. However, the training information will be lost when the Mathematica kernel is closed, so the trained neural net needs to be exported. Exporting allows the trained neural net function to be imported into this or other files. The file can be exported in either wmlf or wlnet formats, both demonstrated below.

VIII. Deployment to Excel

Function-Based Model with Mathematica Link for Excel

Mathematica Link for Excel is a Wolfram application that allows sharing of the Mathematica kernel in Excel. Thus any function, either native to Mathematica or user-defined, can be called from Excel. Deployment of Mathematica functions to Excel allows functions f2 and f4 to be used in Excel-based data maps to other software. For example, Reference 2 shows that once the functions are active in Excel, they can be called from chemical processing software such as Aspen Plus or CHEMCAD.
The code below shows the procedure for calling the co-current membrane cascade function f4 from Excel. The Link application is started first, followed by a command to share the kernel:
The function f4 produces five outputs, including the permeate to retentate split ratio, the retentate and permeate mole fractions, and the retentate and permeate molar flow rates. While we could call f4 directly in Excel, splitting the output using the Part function [[...]] makes it easier to use the results in excel. The five different outputs are defined as individual functions below.
For testing purposes, the following function calls are made to compare the results to Excel:
A screenshot of the Excel worksheet with calls to Mathematica-linked functions is shown below.

Neural Net Model with Mathematica Link for Excel

Calling the function f4 in Mathematica takes about five seconds on the computer used in this document. Calling the function five times as shown above takes five times longer, or about 25 seconds. A similar time is observed when the functions are called from Excel. This time is significantly reduced by using the trained neural net function. The execution time for the five-function-call cell below is less than 10 msec. The five outputs from f4 are first each defined as separate functions, tested, then called from Excel.
A screenshot of the Excel worksheet with calls to Mathematica-linked neural net functions is shown below.

Neural Net Model with Wolfram Cloud Connector

Wolfram Cloud Connector is another Wolfram application that allows sharing of the Mathematica kernel in Excel. Any function deployed to the Wolfram cloud can be called from Excel and used in Excel-based data maps to other software. The advantage of this method of deployment is that an active Mathematica kernel on the local PC is not required.
The code below shows the procedure for calling the co-current membrane cascade function f4 from Excel. The Wolfram CloudBase is set first, followed by a command to connect and allow access to the specific user base.
As before, there are five outputs and the five functions for each of the outputs were defined in the previous section (splitnet2, xnet2, ynet2, nPnet2, and nRnet2). Those five functions were used here to define the API functions for cloud deployment.
Here are some function tests to make sure the API functions are working correctly before deployment.
The functions are then deployed and called from Excel as shown in the screenshot below.
A screenshot of the Excel worksheet with calls to Cloud-deployed neural net functions is shown below.

IX. References

1. Andrew Biaglow, Sam Cowart, Simuck Yuk, Enoch Nagelli, and Corey James, 2024, “3. Improved Membrane Unit in Mathematica Linked to CHEMCAD”, Mendeley Data, Ver. 2, https://doi.org/10.17632/nz7p8bhhs3.2.
2. S.V. Cowart, S.F. Yuk, C.M. James, E.A. Nagelli, A.I. Biaglow, 2024, “Augmented and Cloud Computing With Chemical Process Simulators,” Software: Practice and Experience, https://doi.org/10.1002/spe.3397
3. J.D. Seader, E.J. Henley, and D.K. Roper, Separation Process Principles, New York: Wiley, 2011, pp. 518-519.

CITE THIS NOTEBOOK

Functional and machine learning models for single stage membrane and co-current membrane cascade​
by Samuel V. Cowart, Simuck F. Yuk, Corey M. James, Enoch A. Nagelli, and Andrew I. Biaglow​
Wolfram Community, STAFF PICKS, January 7, 2025
​https://community.wolfram.com/groups/-/m/t/3352664