#! /bin/sh
######################################
##### Model Transformation Tools #####
######################################
# Bourne shell script: ode2odes_r2c
# Reduce ordinary differential equations to differential-algebraic
# equations solution in the form of a c program.
# Euler integration of the state is included.
# NB Arrays should are defined to be one larger than expected
# - the 0 element is not used.
# Copyright (c) P.J.Gawthrop 1997.
###############################################################
## Version control history
###############################################################
## $Id$
## $Log$
## Revision 1.1 1997/01/21 22:54:54 peterg
## Initial revision
##
###############################################################
# Inform user
echo Creating $1_odes.c
# Remove the old log file
rm -f ode_r2c.log
# Use reduce to accomplish the transformation
reduce >ode2odes_r2c.log << EOF
%Read the reduce definitions file
in "$1_def.r";
%Set up the number of argument variables to zero in case the user has forgotten
MTTNVar := 0;
%Read the symbolic parameters file
in "$1_sympar.r";
ON BigFloat, NumVal;
PRECISION 16; %Compatible with Matlab
%OFF Nat;
ON NERO; % Suppress zero elements
%Generate the Header part
OUT "$1_odes.c";
write "/*"$
write "Program to solve ode for system $1"$
write "NB Arrays are defined to be one larger than expected"$
write " - the 0 element is not used."$
write "File $1_odes.c"$
write "Generated by MTT"$
write "*/"$
write " "$
%Program heading
write "#include <stdio.h>"$
write "#include ""$1_sympar.h"" "$
write "/* Declare standard arrays */"$
write "double y[", MTTNy+1, "]; /* $1_ode output */"$
write "double dx[", MTTNx+1, "]; /* $1_ode state derivative */"$
write "double x[", MTTNx+1, "]; /* $1_ode state */"$
write "double u[", MTTNu+1, "]; /* $1_ode input */"$
%External (global) variable list
write "/* External (global) variable list */ "$
IF MTTNvar>0 THEN
BEGIN
FOR i := 1:MTTNvar DO
IF numberp(MTTVar(i,1))
THEN
BEGIN
% Do nowt
END
ELSE
BEGIN
write "double ", MTTVar(i,1), ";"$
END$
END$
write "/* Files */ "$
write " FILE *fopen(), *fps, *fpso;"$
write "main()"$
write " "$
write "{"$
write "/* Counters etc*/ "$
write " double time;"$
write " double dt;"$
write " int i;"$
write " int k;"$
write "/*functions */ "$
write " extern $1_numpar();"$
%Open the output files
write "/* %Open the output file */"$
write "fps = fopen(""$1_odes.m"", ""w""); "$
write "fpso = fopen(""$1_odeso.m"", ""w""); "$
%Set up user-defined constants
write "/* Set up user-defined constants */"$
write " $1_numpar();"$
%Initialise main (Euler) integration loop
write "/* Initialise main (Euler) integration loop */"$
write " time = 0;"$
write " dt = DT/STEPFACTOR;"$
write " for (i=1; i<=", MTTNx, "; i++)"$
write " x[i] = 0.0;"$
%Set up system inputs
write "/* Set up system inputs */"$
write " for (i=1; i<=", MTTNu, "; i++)"$
write " u[i] = 1.0;"$
write " $1_input(0.0);"$
write " fprintf(fps, ""function data = ", "$1_odes \n"");"$
write " fprintf(fps, ""data = [\n"");"$
write " fprintf(fpso, ""function data = ", "$1_odeso \n"");"$
write " fprintf(fpso, ""data = [\n"");"$
% Compute the first output
write "/* Compute the first output */"$
write " $1_ode(y,dx,x,u);"$
%Main (Euler) integration loop
write "/* Main (Euler) integration loop */"$
write " while (time<LAST)"$
write " {"$
%Write to output to file
write "/* Write to output file */"$
write " fprintf(fpso, ""%5.4f "",time);"$
write " for (i=1; i<=", MTTNy, "; i++)"$
write " fprintf(fpso, ""%5.4f "", y[i]);"$
write " fprintf(fpso, ""\n"");"$
%Write to state to file
write "/* Write to state file */"$
write " fprintf(fps, ""%5.4f "",time);"$
write " for (i=1; i<=", MTTNx, "; i++)"$
write " fprintf(fps, ""%5.4f "", x[i]);"$
write " fprintf(fps, ""\n"");"$
% Inner integration loop
write "/* Inner integration loop */"$
write " for (k=1; k<=STEPFACTOR; k++)"$
write " {"$
%Set up system inputs
write " /* Set up system inputs */"$
write " $1_input(time);"$
write " for (i=1; i<=", MTTNx, "; i++)"$
write " x[i] = x[i] + dx[i]*dt;"$
write " $1_ode();"$
write " time = time + dt;"$
write " }"$
write " }"$
%Write to files
%Write to output to file
write "/* Write to output file */"$
write " fprintf(fpso, ""%5.4f "",time);"$
write " for (i=1; i<=", MTTNy, "; i++)"$
write " fprintf(fpso, ""%5.4f "", y[i]);"$
write " fprintf(fpso, ""\n"");"$
%Write to state to file
write "/* Write to state file */"$
write " fprintf(fps, ""%5.4f "",time);"$
write " for (i=1; i<=", MTTNx, "; i++)"$
write " fprintf(fps, ""%5.4f "", x[i]);"$
write " fprintf(fps, ""\n"");"$
write " fprintf(fps, ""];\n"");"$
write " fprintf(fpso, ""];\n"");"$
write "}"$
SHUT "$1_odes.c";
EOF