Hi, I need to model a countercurrent shell-tube exchanger with slightly more accuracy than usual. The model has to respond to flows and feed temperature changes on either input, and has to equate to an LMTD model under static conditions. The dynamics will come from the thermal mass of the exchanger elements and contents, and the transport behaviour inside it.

My default will be to do it as a finite element model, I'm just wondering whether anyone has a reference to a 'clever' alternative approach. TIA

That depends on the tools you have available to you. If you have a good finite element package, then you will have a fairly easy time just plugging it in to that. If you are going to code it all on your own, you might consider just using a finite difference model, basically taking the differential equations and solving using perhaps a runge kutta method.

You could also work out the full analytical solution. It would take a bit more work up front, but would require less computer time later.

A long time ago when I was in power generation I saw how the LMTD equations were developed. It basically divided the length of the counter flow heat exchanger into smaller and smaller sections. I don't see how finite element analysis will help.

Wouldn't a non-linear differential equation do the job? Then use system identification to find the coefficients of the non-linear differential equation? This would required the control and response data as a function of time.

Usually I'm only interested in one or two inputs and outputs, the rest of the responses don't matter. In that case, you can generally just use a couple of deadtime/lag models and that's close enough. Here I need to have correct responses (dynamic and steady state) for several interacting variables. The device is actually a gas/solid contactor/HE - but that's another story.

A long time ago when I was in power generation I saw how the LMTD equations were developed. It basically divided the length of the counter flow heat exchanger into smaller and smaller sections. I don't see how finite element analysis will help.

Wouldn't a non-linear differential equation do the job? Then use system identification to find the coefficients of the non-linear differential equation? This would required the control and response data as a function of time.

Yes, a PDE might be the way to go. I've spent so many years getting messrs Runge & Kutta to do the work for me that it's going to take some effort to get back to that level. :-) Thanks for the input guys.

This is not quite what you are looking for but it is a start. I wrote this to find the coefficients that minimize the minimum integrated absolute error. You would need to modify this

// Find the coefficients that minimize the minimum integrate absolute error

T=3D0.01; // SAMPLE INTERVAL tLength=3D60; // LENGTH IN TIME UNITS N=3Dround(tLength/T); // NUMBER OF INTERVALS r=3D1; // REFERENCE OR SET POINT

First r can no longer be a constant. You would need to replace r by CO(t-DeadTime) in the differential equation. In this case it is the forcing function. DeadTime will be one of those variables optimized by the optim function. The r used in calculating the iae needs to be replaced with PV(t). PV(t) and CO(t) are recorded data like from a RSTrend or something similar.

function ydot=3DDifEQ3(t,y,A) // SOLVE FOR THIRD ORDER ydot=3D[y(2);y(3);-A(2)*y(3)-A(1)*y(2)-y(1)+r]; // DIFFERENTIAL EQUATION endfunction

The DifEQ3 function my need to be expanded to take into account all the different states you think you need. In this case there are three simple first order differential equations but they can be non- linear. As I said above the r need to be replaced by CO(t- DeadTime). If you have two inputs then you will need yet another function and dead time.

The evaluation function minimizes the IAE. The optim will work for your purposes but usually the Levenberg-Marquardt method is used for system identification. Scilab's lsqrsolve function does this. The lsqrsolve uses the array of errors(t)=3DPV(t)-EV(t). EV(t) is the estimated process value returned by the differential equation.

This part would need to be replaced with a call to lsqrsolve which minimizes the sum of squared errors. It would also need to read the PV and CO trend data first.

I have a program that can ID FOPDT and SOPDT written in Scilab. It takes about 30 seconds to run. It also compares the results of the FOPDT and SOPDT and chose the model with the lowest error and calculates the gain. It took about 3 pages to handle all the details. That isn't bad.

I see only two big hurdles to over come.

Determining the general model. I think a FOPDT or SOPDT will work but the gain is non-linear and will vary with temperature. There will need to be a coefficient(s) that will take into account how the gain changes with temperature.

Exciting the system to get good PV(t) and CO(t). Realistically it can't be done with just a simple step in the CO. The CO must reach a steady state at least to levels and maybe three or more in your cases. This is necessary to calculate the gain and offset ( ambient temperature ) There must also be plenty of transitions. Transitions are good for finding time constants and dead times. If you have two inputs then changing the both of them in a manner that will yield useful data will be tricky. This will probably require changing one input at a time.

If you don't want to mess with this you can probably buy an off the shelf auto tuning program but I don't know how well they handle non- linear gains or the general form of the model or two inputs.

This is not quite what you are looking for but it is a start. I wrote this to find the coefficients that minimize the minimum integrated absolute error. You would need to modify this

// Find the coefficients that minimize the minimum integrate absolute error

T=0.01; // SAMPLE INTERVAL tLength=60; // LENGTH IN TIME UNITS N=round(tLength/T); // NUMBER OF INTERVALS r=1; // REFERENCE OR SET POINT

First r can no longer be a constant. You would need to replace r by CO(t-DeadTime) in the differential equation. In this case it is the forcing function. DeadTime will be one of those variables optimized by the optim function. The r used in calculating the iae needs to be replaced with PV(t). PV(t) and CO(t) are recorded data like from a RSTrend or something similar.

function ydot=DifEQ3(t,y,A) // SOLVE FOR THIRD ORDER ydot=[y(2);y(3);-A(2)*y(3)-A(1)*y(2)-y(1)+r]; // DIFFERENTIAL EQUATION endfunction

The DifEQ3 function my need to be expanded to take into account all the different states you think you need. In this case there are three simple first order differential equations but they can be non- linear. As I said above the r need to be replaced by CO(t- DeadTime). If you have two inputs then you will need yet another function and dead time.

The evaluation function minimizes the IAE. The optim will work for your purposes but usually the Levenberg-Marquardt method is used for system identification. Scilab's lsqrsolve function does this. The lsqrsolve uses the array of errors(t)=PV(t)-EV(t). EV(t) is the estimated process value returned by the differential equation.

This part would need to be replaced with a call to lsqrsolve which minimizes the sum of squared errors. It would also need to read the PV and CO trend data first.

I have a program that can ID FOPDT and SOPDT written in Scilab. It takes about 30 seconds to run. It also compares the results of the FOPDT and SOPDT and chose the model with the lowest error and calculates the gain. It took about 3 pages to handle all the details. That isn't bad.

I see only two big hurdles to over come.

Determining the general model. I think a FOPDT or SOPDT will work but the gain is non-linear and will vary with temperature. There will need to be a coefficient(s) that will take into account how the gain changes with temperature.

Exciting the system to get good PV(t) and CO(t). Realistically it can't be done with just a simple step in the CO. The CO must reach a steady state at least to levels and maybe three or more in your cases. This is necessary to calculate the gain and offset ( ambient temperature ) There must also be plenty of transitions. Transitions are good for finding time constants and dead times. If you have two inputs then changing the both of them in a manner that will yield useful data will be tricky. This will probably require changing one input at a time.

If you don't want to mess with this you can probably buy an off the shelf auto tuning program but I don't know how well they handle non- linear gains or the general form of the model or two inputs.

Peter Nachtwey

Peter, Thanks for taking the time. I'll read through this and post some comments in a few days.

PolyTech Forum website is not affiliated with any of the manufacturers or service providers discussed here.
All logos and trade names are the property of their respective owners.