Home > m > timetools > statespacefunctions > ltpda_ss_modify_old.m

ltpda_ss_modify_old

PURPOSE ^

OLD VERSION

SYNOPSIS ^

function varargout = ltpda_ss_modify(varargin)

DESCRIPTION ^

 OLD VERSION

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 DESCRIPTION: ltpda_ss_modify enables to modify parameters and retrieve
 system data in numrical or parametric form, eventually take the
 jacobian or hessian of its fields.

 CALL: varargout = ltpda_ss_modify(varargin)
 [subystem] = paramSubstitute(subystem, operations)
 [subystem] = paramSubstitute(subystem, operations)
 [Subystem_jac, subystem] = paramSubstitute(subystem, operations)
 [Subystem_hess, Subystem_jac, subystem] = paramSubstitute(subystem, operations)
 
 INPUTS: subystem - plist describing a subsystem (refer to ltpda_ss_check for details)
 operations - plist of elements for numerical values the of parameters

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function varargout = ltpda_ss_modify(varargin)
0002 % OLD VERSION
0003 %
0004 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0005 %
0006 % DESCRIPTION: ltpda_ss_modify enables to modify parameters and retrieve
0007 % system data in numrical or parametric form, eventually take the
0008 % jacobian or hessian of its fields.
0009 %
0010 % CALL: varargout = ltpda_ss_modify(varargin)
0011 % [subystem] = paramSubstitute(subystem, operations)
0012 % [subystem] = paramSubstitute(subystem, operations)
0013 % [Subystem_jac, subystem] = paramSubstitute(subystem, operations)
0014 % [Subystem_hess, Subystem_jac, subystem] = paramSubstitute(subystem, operations)
0015 %
0016 % INPUTS: subystem - plist describing a subsystem (refer to ltpda_ss_check for details)
0017 % operations - plist of elements for numerical values the of parameters
0018 
0019 % 'NOT_USED_PARAM_NAMES' cell array of parameter discarted in this
0020 % analysis, substituted by their standard numerical value
0021 % 'NOT_USED_PARAM_ORDER' array of 1 / 0 depending whether parameter are
0022 % effectively discarted or not
0023 %
0024 % 'SET_PARAM_NAMES' Cell array of parameter names for numerical substitutions with user defined values
0025 % 'SET_PARAM_VALUES' Array of parameter means for numerical substitutions
0026 % 'SET_PARAM_SIGMAS' Array of parameter standard deviation for numerical substitutions (unit may be none!!!)
0027 % 'SET_PARAM_VALUES_ORDER' Array of 1 / 0 depending whether parameter means are effectively substituted or not
0028 % 'SET_PARAM_SIGMAS_ORDER' Array of 1 / 0 depending whether parameter means are effectively substituted or not
0029 %
0030 % 'JAC_PARAM_NAMES' same for parameters used for evaluating a jacobian
0031 % 'JAC_PARAM_ORDER'
0032 % 'HESS_PARAM_NAMES' same for parameters used for evaluating a hessian
0033 % 'HESS_PARAM_ORDER'
0034 % 'SYMBOLIC_OUTPUT_ALLOWED' if 0 checks output tensors are double and not symbollic arrays
0035 %
0036 % OUTPUTS: subystem - plist describing the subsystem (refer to test_subsys
0037 % for details) with substitutions / modifications done
0038 % System_jac - plist of plists (subsystems) giving the jacobian of the
0039 % matrices in regard with the differenciation parameters, with the
0040 % specified substitutions made
0041 % System_hess - plist of plists (subsystems) giving the hessian of the
0042 % matrices in regard with the differenciation parameters, with the
0043 % specified substitutions made
0044 %
0045 % in both cases differenciated array are the same size for the first two
0046 % dimensions, the third 8resp. 4rth) one being for the direction of the
0047 % differenciation
0048 %
0049 % PARAMETERS:
0050 % ================ 'modif' field ================
0051 % A: ...
0052 % ================ 'derivatives' field ================
0053 % ================ 'derivatives2' field ================
0054 %
0055 % VERSION: $Id: ltpda_ss_modify_old.html,v 1.1 2008/03/01 12:29:32 hewitson Exp $
0056 %
0057 % HISTORY: dd-mm-yyyy A Grynagier
0058 % Creation
0059 %
0060 % TO DO :
0061 % display in window parametric substitutions and differenciation proceeded
0062 % Define user scenarii to see ensure output is satisfying
0063 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0064 
0065 %% standard calls for LTPDA function data
0066 
0067 ALGONAME = mfilename;
0068 VERSION = '$Id: ltpda_ss_modify_old.html,v 1.1 2008/03/01 12:29:32 hewitson Exp $';
0069 display(['starting ' mfilename]);
0070 
0071 if isequal( varargin{1}, 'Version')
0072     listMeta = VERSION;
0073     return;
0074 elseif isequal(varargin{1}, 'Params')
0075     listMeta = plist();
0076     return;
0077 elseif isequal(varargin{1}, 'Category')
0078     listMeta = 'To be defined / State Space / Time domain analysis';
0079     return;
0080 end
0081 
0082 %% Retriveing option fields for operations
0083 
0084 if length(varargin)==2
0085     System = varargin{1};
0086     modif = varargin{2};
0087     %add loop in case modif is an array of plists
0088 else
0089     error('wrong number of inputs!')
0090 end
0091 % fields from user options
0092 NotUsedParamNamesIni = find(modif,'NOT_USED_PARAM_NAMES') ;
0093 NotUsedParamOrderIni = find(modif,'NOT_USED_PARAM_ORDER') ;
0094 SetParamNamesIni = find(modif,'SET_PARAM_NAMES') ;
0095 SetParamValuesIni = find(modif,'SET_PARAM_VALUES') ;
0096 SetParamSigmasIni = find(modif,'SET_PARAM_SIGMAS') ;
0097 SetParamValuesOrderIni = find(modif,'SET_PARAM_VALUES_ORDER') ;
0098 SetParamSigmasOrderIni = find(modif,'SET_PARAM_SIGMAS_ORDER') ;
0099 JacParamNamesIni = find(modif,'JAC_PARAM_NAMES') ;
0100 JacParamOrderIni = find(modif,'JAC_PARAM_ORDER') ;
0101 HessParamNamesIni = find(modif,'HESS_PARAM_NAMES') ;
0102 HessParamOrderIni = find(modif,'HESS_PARAM_ORDER') ;
0103 SymbolicOutputAllowed = find(modif,'SYMBOLIC_OUTPUT_ALLOWED') ;
0104 
0105 % subsystem data matrices
0106 SystemName = find(System,'NAME');
0107 % SystemTimestep = find(System,'TIMESTEP');
0108 % SystemXini = find(System,'XINI');
0109 SystemAMat = find(System,'AMAT');
0110 SystemBMats = find(System,'BMATS');
0111 SystemCMat = find(System,'CMAT');
0112 SystemDMats = find(System,'DMATS');
0113 % subsystem parametric data
0114 SystemParamNames = find(System,'PARAMNAMES');
0115 SystemParamValues = find(System,'PARAMVALUES');
0116 SystemParamSigmas = find(System,'PARAMSIGMAS');
0117 SystemNbParams = length(SystemParamNames);
0118 
0119 %% building updated lists
0120 % Not Used
0121 if isequal(NotUsedParamOrder,'ALL')
0122     NotUsedParamNames = NotUsedParamNamesIni;
0123 elseif isequal(NotUsedParamOrder,'NO')
0124     NotUsedParamNames = {};
0125 else
0126     j=1;
0127     for i=1:legnth(NotUsedParamNamesIni)
0128         if NotUsedParamNamesIni(i)
0129             NotUsedParamNames{j} = NotUsedParamNamesIni{i};
0130             j = j+1;
0131         end
0132     end
0133 end
0134 % Set
0135 if isequal(SetParamOrder,'ALL')
0136     SetParamNames = SetParamNamesIni;
0137 elseif isequal(SetParamOrder,'NO')
0138     SetParamNames = {};
0139 else
0140     j=1;
0141     for i=1:legnth(SetParamNamesIni)
0142         if SetParamNamesIni(i)
0143             SetParamNames{j} = SetParamNamesIni{i};
0144             j = j+1;
0145         end
0146     end
0147 end
0148 % Jac
0149 if isequal(JacParamOrder,'ALL')
0150     JacParamNames = JacParamNamesIni;
0151 elseif isequal(JacParamOrder,'NO')
0152     JacParamNames = {};
0153 else
0154     j=1;
0155     for i=1:legnth(JacParamNamesIni)
0156         if JacParamNamesIni(i)
0157             JacParamNames{j} = JacParamNamesIni{i};
0158             j = j+1;
0159         end
0160     end
0161 end
0162 % Hess
0163 if isequal(HessParamOrder,'ALL')
0164     HessParamNames = HessParamNamesIni;
0165 elseif isequal(HessParamOrder,'NO')
0166     HessParamNames = {};
0167 else
0168     j=1;
0169     for i=1:legnth(HessParamNamesIni)
0170         if HessParamNamesIni(i)
0171             HessParamNames{j} = HessParamNamesIni{i};
0172             j = j+1;
0173         end
0174     end
0175 end
0176 nNotUsed = length(NotUsedParamNames);
0177 nSet = length(SetParamNames);
0178 nJac = length(JacParamNames);
0179 nHess = length(HessParamNames);
0180 
0181 %% building param tables
0182 InSet = zeros(SystemNbParams,1);
0183 InJac = zeros(SystemNbParams,1);
0184 InHess = zeros(SystemNbParams,1);
0185 InNotUsed = zeros(SystemNbParams,1);
0186 for i=1:SystemNbParams
0187     for j=1:nNotUsed
0188         if isequal(SystemParamNames{i},NotUsedParamNames{j})
0189             InNotUsed(i) = j;
0190         end
0191         if isequal(SystemParamNames{i},SetParamNames{j})
0192             InSet(i) = j;
0193         end
0194         if isequal(SystemParamNames{i},JacParamNames{j})
0195             InJac(i) = j;
0196         end
0197         if isequal(SystemParamNames{i},HessParamNames{j})
0198             InHess(i) = j;
0199         end
0200     end
0201 end
0202 
0203 %% looking for non exisitng but modified parameter names
0204 if not(length(NotUsedParamNames)==sum(InNotUsed>0))
0205     msg = 'error with notusedparameter list because one does not exist in system.';
0206     error(msg);
0207 elseif not(length(SetParamNames)==sum(InSet>0))
0208     msg = 'error with setparameter list because one does not exist in system.';
0209     error(msg);
0210 elseif not(length(JacParamNames)==sum(InJac>0))
0211     msg = 'error with jacparameter list because one does not exist in system.';
0212     error(msg);
0213 elseif not(length(HessParamNames)==sum(InHess>0))
0214     msg = 'error with hessparameter list because one does not exist in system.';
0215     error(msg);
0216 end
0217 
0218 %% looking of unused but used parameters
0219 for i=1:nNotUsed
0220     for j =1:nSet
0221         if isequal(NotUsedParamNames{i},SetParamNames{j})
0222             msg = 'error because parameter set was supposed to be unused';
0223             error(msg);
0224         end
0225     end
0226     for j =1:nJac
0227         if isequal(NotUsedParamNames{i},JacParamNames{j})
0228             msg = 'error because parameter for jacobian was supposed to be unused';
0229             error(msg);
0230         end
0231     end
0232     for j =1:nHess
0233         if isequal(NotUsedParamNames{i},HessParamNames{j})
0234             msg = 'error because parameter for hessian was supposed to be unused';
0235             error(msg);
0236         end
0237     end
0238 end
0239 
0240 %% Substituting unused parameters
0241 for i=1:nNotUsed
0242     SystemAMat{1} = subs(SystemAMat{1}, NotUsedParamNames(i) , NotUsedParamValues(i));
0243     SystemCMat{1} = subs(SystemCMat{1}, NotUsedParamNames(i) , NotUsedParamValues(i));
0244     for j=1:length(SystemBMats)
0245         SystemBMats{j} = subs(SystemBMats{j}, NotUsedParamNames(i) , NotUsedParamValues(i));
0246         SystemDMats{j} = subs(SystemDMats{j}, NotUsedParamNames(i) , NotUsedParamValues(i));
0247     end
0248 end
0249 
0250 %% modifying parameter tables in System (setparam list)
0251 
0252 if not(isequal(SetParamValuesOrder,'NO')) && not(isequal(SetParamValuesOrder,[]))
0253     % check all parameters to be set exist in the initial System
0254     k = 1;
0255     for i=1:length(SetParamNames)
0256         found = 0;
0257         for j=1:length(SetParamNames)
0258             if isequal(SetParamNames{i},SystemParamNames{j})
0259                 found = 1;
0260                 SetParamNames2{k} =  SetParamNames{j};
0261                 SetParamValues2(k) =  SetParamValues(j);
0262                 SetParamSigmas2(k) =  SetParamSigmas(j);
0263                 k = k+1;
0264             end
0265         end
0266         if found == 0 && (SetParamValuesOrder(i) || SetParamSigmasOrder(i) || isequal(SetParamValuesOrder,'ALL') || isequal(SetParamSigmasOrder,'ALL'))
0267             msg = ['error because parameter ', SetParamNames{i}, ' was not found in the system ',SystemName];
0268             error(msg)
0269         end
0270     end
0271     %proceeding substitutions into standard values - not into symbollic
0272     %expressions yet
0273     for i=1:length(SetParamNames)
0274         for j=1:length(SystemParamNames)
0275             if isequal(SetParamNames{i},SystemParamNames{j})
0276                 if SetParamValuesOrder(i)||isequal(SetParamValuesOrder,'ALL')
0277                     SystemParamValues(i) = SetParamValues(j);
0278                 end
0279                 if SetParamSigmasOrder(i)||isequal(SetParamSigmasOrder,'ALL')
0280                     SystemParamSigmas(i) = SetParamSigmas(j);                    
0281                 end
0282             end
0283         end
0284     end
0285 end
0286 
0287 %% building jacobians - symbollic expression
0288 % building vector for differenciation
0289 JacParamIndices = zeros(length(JacParamNames));
0290 for i=1:length(JacParamNames)
0291     for j=1:length(SetParamNames)
0292         if isequal(JacParamNames{i},SetParamNames{j})&&JacParamIndices(i)==0
0293             JacParamIndices(i)=j;
0294         end
0295     end
0296 end
0297 if not(isequal(JacParamOrder, 'NO'))
0298     i = 1;
0299     JacParamNames2 = {};
0300     for j=1:length(JacParamNames)
0301         if isequal(JacParamOrder(i),1) || isequal(JacParamOrder(i),'ALL')
0302             JacParamNames2{i} =  JacParamNames{j};
0303             JacParamValues2(i) =  SetParamValues(JacParamIndices(j));
0304             JacParamSigmas2(i) =  SetParamSigmas(JacParamIndices(j));
0305             i = i+1;
0306         end
0307     end
0308     if  isequal(JacParamNames2, {}) 
0309         msg = 'warning the parameter list for differenciation is empty';
0310         display(msg);
0311     end
0312     % creating symbollic variables vector
0313     JacParamVar = zeros(length(JacParamNames2),1);
0314     for j = 1:length(JacParamNames2)
0315         nameVar = eval(['JacParamNames2{', num2str(j),'}']);
0316 %         cmd = [nameVar ,'=sym(', nameVar ,');']; Declaration not needed !!
0317 %         eval(cmd);
0318         cmd2 = ['JacParamVar(',num2str(j),')=',nameVar];
0319         eval(cmd2);
0320     end
0321     % differenciating
0322     JacSystemAMat{1} = ltpda_jacobian(SystemAMat{1}, JacParamVar,1);
0323     JacSystemCMat{1} = ltpda_jacobian(SystemCMat{1}, JacParamVar,1);
0324     for j=1:length(SystemBMats)
0325         JacSystemBMats{j} = ltpda_jacobian(SystemBMats{j}, JacParamVar,1);
0326         JacSystemDMats{j} = ltpda_jacobian(SystemDMats{j}, JacParamVar,1);
0327     end
0328 end
0329 
0330 %% building hessians - symbollic expression
0331 HessParamIndices = zeros(length(HessParamNames));
0332 for i=1:length(HessParamNames)
0333     for j=1:length(SetParamNames)
0334         if isequal(HessParamNames{i},SetParamNames{j})&&HessParamIndices(i)==0
0335             HessParamIndices(i)=j;
0336         end
0337     end
0338 end
0339 if not(isequal(HessParamOrder, 'NO'))
0340     i = 1;
0341     HessParamNames2 = {};
0342     for j=1:length(HessParamNames)
0343         if isequal(HessParamOrder(i),1) || isequal(HessParamOrder(i),'ALL')
0344             HessParamNames2{i} =  HessParamNames{j};
0345             HessParamValues2(i) =  SetParamValues(HessParamIndices(j));
0346             HessParamSigmas2(i) =  SetParamSigmas(HessParamIndices(j));
0347             i = i+1;
0348         end
0349     end
0350     if  isequal(HessParamNames2, {}) 
0351         msg = 'warning the parameter list for second order differenciation is empty';
0352         display(msg);
0353     end
0354     % creating symbollic variables vector
0355     HessParamVar = zeros(length(HessParamNames2),1);
0356     for j = 1:length(HessParamNames2)
0357         nameVar = eval(['HessParamNames2{', num2str(j),'}']);
0358 %         cmd = [nameVar ,'=sym(', nameVar ,');']; Declaration not needed !!
0359 %         eval(cmd);
0360         cmd2 = ['HessParamVar(',num2str(j),')=',nameVar];
0361         eval(cmd2);
0362     end
0363 
0364     % differenciating
0365     HessSystemAMat{1} = ltpda_jacobian(SystemAMat{1}, HessParamVar,2);
0366     HessSystemCMat{1} = ltpda_jacobian(SystemCMat{1}, HessParamVar,2);
0367     for j=1:length(SystemBMats)
0368         HessSystemBMats{j} = ltpda_jacobian(SystemBMats{j}, HessParamVar,2);
0369         HessSystemDMats{j} = ltpda_jacobian(SystemDMats{j}, HessParamVar,2);
0370     end
0371 end
0372 
0373 %% Numerical value for jacobians
0374 
0375 if not(isequal(JacParamOrder, 'NO')) && not(isequal(JacParamOrder,[]))
0376     JacSystemAMat{1} = subs(SystemAMat{1}, JacParamNames2 , JacParamValues2);
0377     JacSystemCMat{1} = subs(SystemCMat{1}, JacParamNames2 , JacParamValues2);
0378     for j=1:length(SystemBMats)
0379         JacSystemBMats{j} = subs(SystemBMats{j}, JacParamNames2 , JacParamValues2);
0380         JacSystemDMats{j} = subs(SystemDMats{j}, JacParamNames2 , JacParamValues2);
0381     end
0382 end
0383 
0384 %% Numerical value for hessians
0385 
0386 if not(isequal(HessParamOrder, 'NO')) && not(isequal(HessParamOrder,[]))
0387     HessSystemAMat{1} = subs(SystemAMat{1}, HessParamNames2 , HessParamValues2);
0388     HessSystemCMat{1} = subs(SystemCMat{1}, HessParamNames2 , HessParamValues2);
0389     for j=1:length(SystemBMats)
0390         HessSystemBMats{j} = subs(SystemBMats{j}, HessParamNames2 , HessParamValues2);
0391         HessSystemDMats{j} = subs(SystemDMats{j}, HessParamNames2 , HessParamValues2);
0392     end
0393 end
0394 
0395 %% Numerical value for system
0396 
0397 if not(isequal(SetParamOrder, 'NO')) && not(isequal(SetParamValuesOrder,[]))
0398     SystemAMat{1} = subs(SystemAMat{1}, SetParamNames2 , SetParamValues2);
0399     SystemCMat{1} = subs(SystemCMat{1}, SetParamNames2 , SetParamValues2);
0400     for j=1:length(SystemBMats)
0401         SystemBMats{j} = subs(SystemBMats{j}, SetParamNames2 , SetParamValues2);
0402         SystemDMats{j} = subs(SystemDMats{j}, SetParamNames2 , SetParamValues2);
0403     end
0404 end
0405 
0406 if not(SYMBOLIC_OUTPUT_ALLOWED)
0407     display('Forbidden symbollic output case is not handled in the current version' );
0408 end
0409 end

Generated on Fri 22-Feb-2008 23:32:26 by m2html © 2003