1/* 2This subroutine computes state variable equations 3Copyright (C) 1999 Dan Stanger 4 5This library is free software; you can redistribute it and/or modify it 6under the terms of the GNU Library General Public License as published 7by the Free Software Foundation; either version 2 of the License, or (at 8your option) any later version. 9 10This library is distributed in the hope that it will be useful, but 11WITHOUT ANY WARRANTY; without even the implied warranty of 12MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 13Library General Public License for more details. 14 15You should have received a copy of the GNU Library General Public 16License along with this library; if not, write to the Free Software 17Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. 18 19Dan Stanger dan.stanger@ieee.org 20*/ 21/* load in readfile and tree */ 22/* get the filename from the user */ 23/* this code works around problems caused by: 24 matrix.scalar does not simplify properly 25 1 dimensional matrixes do not simplify properly 26 empty matrixes do not remember how they were created 27*/ 28f(m,n,p,q,r,s, MN, MB):=block([l,a_T, a_L,b_T,b_L,c_T,c_L], 29 /* reorder the array a_pt and turn it into the a_T and a_L matrixes */ 30 l:sortem(MN,MB),a_T:first(l),a_L:second(l), 31 b_T:-transpose(a_L).invert(transpose(a_T)), 32 b_L:ident(MN), 33 c_T:ident(MN), 34 c_L:-transpose(b_T), 35 block([ 36 em:matrix(), 37 c_LLC, c_LRC, c_LGC, c_LLR, c_LRR, c_LGR, c_LLG, c_LRG, c_LGG, 38 b_TCL, b_TRL, b_TGL, b_TCR, b_TRR, b_TGR, b_TCG, b_TRG, b_TGG, 39 x:transpose(append( /* state variable vector */ 40 makelist(concat('v,pt_name[i]),i,1,m), 41 makelist(concat('i,pt_name[i]),i,m+n+p+1,m+n+p+q))), 42 u:transpose(append(makelist( 43 if pt_V[i] = false then pt_name[i] else pt_V[i], 44 i,m+n+1,m+n+p), 45 makelist( 46 if pt_V[i] = false then pt_name[i] else pt_V[i], 47 i,m+n+p+q+r+1,m+n+p+q+r+s))), 48 r_L,g_T,se,sx,su,e1,e2,e3,e4], 49 if length(x) = 1 then x:x[1], 50 if length(u) = 1 then u:u[1], 51 submat_m(c_LLC, c_L, 1 .. m, 1 .. q), 52 submat_m(c_LRC, c_L, 1 .. m, (q+1) .. (q+r)), 53 submat_m(c_LGC, c_L, 1 .. m, (q+r+1) .. (q+r+s)), 54 submat_m(c_LLR, c_L, (m+1) .. (m+n), 1 .. q), 55 submat_m(c_LRR, c_L, (m+1) .. (m+n), (q+1) .. (q+r)), 56 submat_m(c_LGR, c_L, (m+1) .. (m+n), (q+r+1) .. (q+r+s)), 57 submat_m(c_LLG, c_L, (m+n+1) .. (m+n+p), 1 .. q), 58 submat_m(c_LRG, c_L, (m+n+1) .. (m+n+p), (q+1) .. (q+r)), 59 submat_m(c_LGG, c_L, (m+n+1) .. (m+n+p), (q+r+1) .. (q+r+s)), 60 61 submat_m(b_TCL, b_T, 1 .. q, 1 .. m), 62 submat_m(b_TRL, b_T, 1 .. q, (m+1) .. (m+n)), 63 submat_m(b_TGL, b_T, 1 .. q, (m+n+1) .. (m+n+p)), 64 submat_m(b_TCR, b_T, (q+1) .. (q+r), 1 .. m), 65 submat_m(b_TRR, b_T, (q+1) .. (q+r), (m+1) .. (m+n)), 66 submat_m(b_TGR, b_T, (q+1) .. (q+r), (m+n+1) .. (m+n+p)), 67 submat_m(b_TCG, b_T, (q+r+1) .. (q+r+s), 1 .. m), 68 submat_m(b_TRG, b_T, (q+r+1) .. (q+r+s), (m+1) .. (m+n)), 69 submat_m(b_TGG, b_T, (q+r+1) .. (q+r+s), (m+n+1) .. (m+n+p)), 70 71 r_L:makelist(pt_name[i],i,m+n+p+q+1,m+n+p+q+r), 72 if length(r_L) > 0 then r_L:apply(diag_matrix,r_L) else r_L:em, 73 g_T:makelist(1/pt_name[i],i,m+1,m+n), 74 if length(g_T) > 0 then g_T:apply(diag_matrix,g_T) else g_T:em, 75 /* compute free response */ 76 e1:mat_unblocker( 77 matrix2([zeromatrix2(mat_nrows_m(c_LRC),mat_ncols_m(b_TRL)),c_LRC], 78 [b_TRL,zeromatrix2(mat_nrows_m(b_TRL),mat_ncols_m(c_LRC))])), 79 e2:block([e:matrix2([b_TRR,r_L], [g_T, c_LRR])], 80 if matrixp(e) then invert(e) else 1/e), 81 e3:mat_unblocker( 82 matrix2( [b_TCR, zeromatrix2(mat_nrows_m(b_TCR),mat_ncols_m(c_LLR))], 83 [zeromatrix2(mat_nrows_m(c_LLR),mat_ncols_m(b_TCR)),c_LLR])), 84 e4:mat_unblocker( 85 matrix2([zeromatrix2(mat_nrows_m(c_LLC),mat_ncols_m(b_TCL)),c_LLC], 86 [b_TCL,zeromatrix2(mat_nrows_m(b_TCL),mat_ncols_m(c_LLC))])), 87 sx:(e1.e2.e3), 88 if e4 # em then sx:sx-e4, /* check no link inductors */ 89 sx:if matrixp(x) then sx.x else x*sx, /* if x is scalar use * mult */ 90 /* compute source response */ 91 e3:mat_unblocker( 92 matrix2( [b_TGR, zeromatrix2(mat_nrows_m(b_TGR),mat_ncols_m(c_LGR))], 93 [zeromatrix2(mat_nrows_m(c_LGR),mat_ncols_m(b_TGR)),c_LGR])), 94 e4:mat_unblocker( 95 matrix2([zeromatrix2(mat_nrows_m(c_LGC),mat_ncols_m(b_TGL)),c_LGC], 96 [b_TGL,zeromatrix2(mat_nrows_m(b_TGL),mat_ncols_m(c_LGC))])), 97 su:(e1.e2.e3), 98 if e4 # em then su:su-e4, 99 su:if matrixp(u) then su.u else u*su, 100 se:sx+su, 101 /* try to remove . if its there */ 102 if not matrixp(se) then subst("*",".",sx+su) else se 103 ) 104)$ 105state():= 106block([filename:read("enter the filename"), se], 107 /* compute the proper tree and use apply to destructure it */ 108 se:apply(f,propertree(filename))) 109$ 110