Home > functions > internal > dsParseModelEquations.m

dsParseModelEquations

PURPOSE ^

PARSEMODELEQUATIONS - parse equations and organize model data in DynaSim model structure

SYNOPSIS ^

function [model,name_map] = dsParseModelEquations(text,varargin)

DESCRIPTION ^

PARSEMODELEQUATIONS - parse equations and organize model data in DynaSim model structure

 Usage:
   model = dsParseModelEquations(STRING,'param',value,...)

 Inputs:
   - STRING (required): one of:
     - string with equations
     - string with name of file containing equations (.eqns or .mech)
   - options (using key/value pairs: 'option1',value1,...):
     - 'namespace': added as prefix to beginning of parameter/etc names (default: '')
     - 'delimiter': separates expressions on same line of model text (default: ';')
   - user-supplied parameter values: ('key',value): name (key) of parameters to
                                     be set and associated user-supplied values
 Outputs:
   - model: DynaSim model structure (see dsCheckModel for details)
   - name_map: useful for namespace-specific substitutions across multiple
     sub-models, see description in dsGenerateModel for more information {name,
     namespace_name, namespace, type}

 Notes:
   - NOTE 1: .eqns files contain fully self contained model equations; .mech
     files define (sub)models that depend on variables linked from elsewhere.
     However, this function does not distinguish between the two.

 Examples:
     model = dsParseModelEquations('dx/dt=3*a*x; x(0)=0','a',0);
     model = dsParseModelEquations('dx/dt=3*a*x, x(0)=0','a',0,'delimiter',',');
     model = dsParseModelEquations('CalciumPump.mech','namespace','HH');
     model = dsParseModelEquations('LIFneuron.eqns');
     model = dsParseModelEquations('a=2; b=2*a; f(x)=b; dx/dt=f(x); x(0)=0; if(x>1)(x=0); current->f(x); monitor f(x); % comments')

   - parsing individual sub-models from specification:
     equations=specification.populations(1).equations;
     [model,map] = dsParseModelEquations(equations,'namespace','pop')
     population_mechanism=specification.populations(1).mechanism_list{1};
     [model,map] = dsParseModelEquations(population_mechanism,'namespace','pop_mech')
     connection_mechanism=specification.connections(1).mechanism_list{1};
     [model,map] = dsParseModelEquations(connection_mechanism,'namespace','pop_pop_mech')

 See also: dsClassifyEquation, dsGenerateModel, dsLocateModelFiles

 Author: Jason Sherfey, PhD <jssherfey@gmail.com>
 Copyright (C) 2016 Jason Sherfey, Boston University, USA

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function [model,name_map] = dsParseModelEquations(text,varargin)
0002 %PARSEMODELEQUATIONS - parse equations and organize model data in DynaSim model structure
0003 %
0004 % Usage:
0005 %   model = dsParseModelEquations(STRING,'param',value,...)
0006 %
0007 % Inputs:
0008 %   - STRING (required): one of:
0009 %     - string with equations
0010 %     - string with name of file containing equations (.eqns or .mech)
0011 %   - options (using key/value pairs: 'option1',value1,...):
0012 %     - 'namespace': added as prefix to beginning of parameter/etc names (default: '')
0013 %     - 'delimiter': separates expressions on same line of model text (default: ';')
0014 %   - user-supplied parameter values: ('key',value): name (key) of parameters to
0015 %                                     be set and associated user-supplied values
0016 % Outputs:
0017 %   - model: DynaSim model structure (see dsCheckModel for details)
0018 %   - name_map: useful for namespace-specific substitutions across multiple
0019 %     sub-models, see description in dsGenerateModel for more information {name,
0020 %     namespace_name, namespace, type}
0021 %
0022 % Notes:
0023 %   - NOTE 1: .eqns files contain fully self contained model equations; .mech
0024 %     files define (sub)models that depend on variables linked from elsewhere.
0025 %     However, this function does not distinguish between the two.
0026 %
0027 % Examples:
0028 %     model = dsParseModelEquations('dx/dt=3*a*x; x(0)=0','a',0);
0029 %     model = dsParseModelEquations('dx/dt=3*a*x, x(0)=0','a',0,'delimiter',',');
0030 %     model = dsParseModelEquations('CalciumPump.mech','namespace','HH');
0031 %     model = dsParseModelEquations('LIFneuron.eqns');
0032 %     model = dsParseModelEquations('a=2; b=2*a; f(x)=b; dx/dt=f(x); x(0)=0; if(x>1)(x=0); current->f(x); monitor f(x); % comments')
0033 %
0034 %   - parsing individual sub-models from specification:
0035 %     equations=specification.populations(1).equations;
0036 %     [model,map] = dsParseModelEquations(equations,'namespace','pop')
0037 %     population_mechanism=specification.populations(1).mechanism_list{1};
0038 %     [model,map] = dsParseModelEquations(population_mechanism,'namespace','pop_mech')
0039 %     connection_mechanism=specification.connections(1).mechanism_list{1};
0040 %     [model,map] = dsParseModelEquations(connection_mechanism,'namespace','pop_pop_mech')
0041 %
0042 % See also: dsClassifyEquation, dsGenerateModel, dsLocateModelFiles
0043 %
0044 % Author: Jason Sherfey, PhD <jssherfey@gmail.com>
0045 % Copyright (C) 2016 Jason Sherfey, Boston University, USA
0046 
0047 %% auto_gen_test_data_flag argin
0048 options = dsCheckOptions(varargin,{'auto_gen_test_data_flag',0,{0,1}},false);
0049 if options.auto_gen_test_data_flag
0050   varargs = varargin;
0051   varargs{find(strcmp(varargs, 'auto_gen_test_data_flag'))+1} = 0;
0052   varargs(end+1:end+2) = {'unit_test_flag',1};
0053   argin = [{text}, varargs]; % specific to this function
0054 end
0055 
0056 model=[];
0057 name_map={};
0058 
0059 % organize optional user-supplied info
0060 % key/value pairs
0061 if nargin>2 % check for at least text input and one key/value pair
0062   keys=varargin(1:2:end); % parameters to set
0063   values=varargin(2:2:end); % values to use
0064 else
0065   keys=[];
0066   values=[];
0067 end
0068 
0069 % set namespace
0070 if ~isempty(keys) && ismember('namespace',keys) % check for user-supplied namespace (i.e., namespace)
0071   namespace=values{ismember(keys,'namespace')}; % user-supplied namespace
0072   if ~isempty(namespace)
0073     namespace=[namespace '_'];
0074   else
0075     namespace='';
0076   end
0077 else
0078   namespace='';
0079 end
0080 
0081 % set delimiter
0082 if ~isempty(keys) && ismember('delimiter',keys) % check for user-supplied delimiter
0083   delimiter = values(ismember(keys,'delimiter')); % user-supplied delimiter
0084 else
0085   delimiter=';';
0086 end
0087 
0088 % error handling for improper input format
0089 if ~ischar(namespace)
0090   error('model "namespace" must be a string.');
0091 end
0092 
0093 if ~ischar(delimiter)
0094   error('expression "delimiter" must be a string.');
0095 end
0096 
0097 %% 1.0 convert text into cell array of strings (one string per line)
0098 % check for DynaSim extensions if input is single \w+ string
0099 
0100 if ischar(text) && isempty(regexp(text,'[^\w.]','once')) && ~any(which(text)) % isempty(regexp(text,'[^\w]','once'))
0101 
0102   %if ischar(text) && ~any(which(text)) && isempty(regexp(text,'[^\w.]','once')) % isempty(regexp(text,'[^\w]','once'))
0103   %if ischar(text) && ~exist(text,'file') && isempty(regexp(text,'[^\w.]','once')) % isempty(regexp(text,'[^\w]','once'))
0104   [~,text]=dsLocateModelFiles(text);
0105   if iscell(text) && ~isempty(text)
0106     text=text{1};
0107   end
0108 end
0109 
0110 % check if input is a filename
0111 if ischar(text) && exist(text,'file')
0112   text = dsReadText(text);
0113 %   [~,name,ext]=fileparts2(text);
0114 %   switch ext
0115 %     case '.m'
0116 %       model=feval(name); % evaluate model-creating function and return model
0117 %       return;
0118 %     case '.mat' % todo: uncomment once dsImportModel supports loading .mat
0119 %       %model=dsImportModel(text);
0120 %       %return;
0121 %   end
0122 %
0123 %   % load equations from file
0124 %   [text,res]=readtext(text,'\n','%'); % text: cell array of strings, one element per line in text file
0125 %
0126 %   % remove all lines without text
0127 %   text=text(res.stringMask);
0128 %
0129 %   % remove leading/trailing white space
0130 %   text=strtrim(text);
0131 %
0132 %   % end each line with semicolon
0133 %   for i=1:length(text)
0134 %     if ~isequal(text{i}(end),';')
0135 %       text{i}(end+1)=';';
0136 %     end
0137 %   end
0138 %
0139 %   % concatenate into a single string
0140 %   text=[text{:}]; % concatenate text from all lines
0141 end
0142 
0143 % split string into cell array of lines delimited by semicolon
0144 if ischar(text)
0145   % remove end-line semicolon if present so split lines are free of all ';'
0146   if text(end)==';'
0147     text=text(1:end-1);
0148   end
0149   
0150   % account for the one exception where ';' does not delimit lines:
0151   % conditional actions with multiple statements (expr1; expr2)
0152   % approach: replace ';' by ',' here then reverse the replacement below
0153   % when storing the action in model.conditionals
0154   pattern='(if\([^;]+\)\s*\([^;\)]+);([^;]+\))'; % if(condiiton)(action1;action2)
0155   replace='$1,$2';
0156   text=regexprep(text,pattern,replace,'ignorecase');
0157   
0158   % now split string into cell array of lines
0159   text = strtrim(regexp(text,delimiter,'split'));
0160 end
0161 if ~iscellstr(text)
0162   error('input not recognized. equations must be provided in single string, cell array of strings, or a text file');
0163 end
0164 
0165 %% 2.0 classify and parse lines; store info in model structure
0166 % use empty struct for Octave compatibility
0167 model.parameters=struct('');
0168 model.fixed_variables=struct('');
0169 model.functions=struct('');
0170 model.monitors=struct('');
0171 model.state_variables={};
0172 model.ODEs=struct('');
0173 model.ICs=struct('');
0174 model.conditionals=struct('');
0175 model.linkers=struct('');
0176 model.comments={};
0177 for index=1:length(text) % loop over lines of text
0178   % organize model data in model structure
0179   line=text{index}; % choose a single expression (lines with multiple expressions have already been split above using regexp-split)
0180   [line,comment]=remove_comment(line); % remove comments
0181   if isempty(line) % e.g., entire line was a comment, or there was nothing there originally
0182     if ~isempty(comment)
0183       model.comments{end+1}=comment;
0184     end
0185     continue;
0186   end
0187   
0188   switch dsClassifyEquation(line,delimiter) % classify
0189     case 'parameter'        % var=(string or number)
0190       rhs=regexp(line,'=(.+)$','tokens','once');
0191       lhs=regexp(line,'^([\w\.]+)\s*=','tokens','once');
0192       lhs{1}=strrep(lhs{1},'.','_'); % e.g., Na.g --> Na_g
0193       name=strtrim(lhs{1}); expression=rhs{1};
0194       model.parameters(1).([namespace name]) = expression;
0195       name_map(end+1,:) = {name,[namespace name],namespace,'parameters'};
0196       if ~isempty(comment)
0197         model.comments{end+1}=sprintf('%s (parameter): %s',[namespace name],comment);
0198       end
0199     case 'fixed_variable'   % var=(expression with grouping or arithmetic), var(#), var([#]), var([# #]), var([#,#]), var(#:#), var(#:end), var([#:#]), var([#:end])
0200       lhs=regexp(line,'^(\w+)\s*=','tokens','once');
0201       if ~isempty(lhs)        
0202         rhs=regexp(line,'=(.+)$','tokens','once');
0203         name=strtrim(lhs{1}); expression=rhs{1};
0204         model.fixed_variables(1).([namespace name]) = expression;
0205         name_map(end+1,:) = {name,[namespace name],namespace,'fixed_variables'};
0206       else
0207         % check for update to fixed variable that is already defined
0208         lhs=regexp(line,'^(.+)\(.*\)\s*=','tokens','once');
0209         name=strtrim(lhs{1});
0210         if isfield(model.fixed_variables(1),[namespace name])
0211           % add update to fixed variable definition
0212           expression=[model.fixed_variables(1).([namespace name]) ';' line];
0213           model.fixed_variables(1).([namespace name]) = expression;
0214         else
0215           warning('failed to set fixed variable.');
0216         end
0217       end
0218       if ~isempty(comment)
0219         model.comments{end+1}=sprintf('%s (fixed_variable): %s',[namespace name],comment);
0220       end
0221     case 'function'         % f(vars)=exression
0222 %       if any(line=='@')
0223 %         line=strrep(line,'@','');
0224 %         %error('model specification error: delete the ''@'' character from all function definitions and try again.');
0225 %       end
0226       name=regexp(line,'^(.+)\(.*\)\s*=','tokens','once');
0227       vars=regexp(line,'\((.+)\)\s*=','tokens','once');
0228       rhs=regexp(line,'=(.+)$','tokens','once');
0229       name=strtrim(name{1});
0230       expression=sprintf('@(%s)%s',vars{1},rhs{1});
0231       model.functions(1).([namespace name]) = expression;
0232       name_map(end+1,:) = {name,[namespace name],namespace,'functions'};
0233       if ~isempty(comment)
0234         model.comments{end+1}=sprintf('%s (function): %s',[namespace name],comment);
0235       end
0236     case 'ODE'              % x'=expression or dx/dt=expression
0237       var=regexp(line,'^d(\w+)/dt\s*=','tokens','once'); % x from dx/dt=
0238       if isempty(var)
0239         var=regexp(line,'^(\w+)''\s*=','tokens','once'); % x from x'=
0240       end
0241       rhs=regexp(line,'=(.+)$','tokens','once');
0242       state_variable=strtrim(var{1}); expression=rhs{1};
0243       model.ODEs(1).([namespace state_variable])=expression;
0244       if ~ismember([namespace state_variable],model.state_variables)
0245         name_map(end+1,:) = {state_variable,[namespace state_variable],namespace,'state_variables'};
0246         model.state_variables{end+1}=[namespace state_variable];
0247       end
0248       if ~isempty(comment)
0249         model.comments{end+1}=sprintf('d/dt %s (ODE): %s',[namespace state_variable],comment);
0250       end
0251     case 'IC'               % x(0)=expression
0252       var=regexp(line,'^(\w+)\(','tokens','once');
0253       rhs=regexp(line,'=(.+)$','tokens','once');
0254       state_variable=strtrim(var{1}); expression=rhs{1};
0255       model.ICs(1).([namespace state_variable])=expression;
0256       if ~isempty(comment)
0257         model.comments{end+1}=sprintf('%s(0) (IC): %s',[namespace state_variable],comment);
0258       end
0259     case 'monitor'          % monitor f=(expression or function)
0260       % split list of monitors
0261       lines=strtrim(regexp(line,',','split'));
0262       % loop over monitors in list
0263       for l=1:length(lines)
0264         % process this monitor
0265         line=lines{l};
0266         
0267         % split left and right parts of monitor
0268         lhs=regexp(line,'^monitor ([\w,@\s\.]+)','tokens','once');
0269         if isempty(lhs)
0270           lhs=regexp(line,'([\w,@\s\.]+)','tokens','once');
0271         end
0272         rhs=regexp(line,'=(.+)$','tokens','once');
0273         
0274         % expand list of monitor names (e.g., monitor iNa.I, iK.I)
0275         names=strtrim(regexp(lhs{1},',','split'));
0276         for i=1:length(names) % loop over list of monitors on this line
0277           name=names{i};
0278           % process special monitors (those including '.', e.g., v.spikes(0))
0279           % todo: clean up or generalize this procedure...
0280           if any(name=='.')
0281             % check for numeric monitor argument
0282             arg=regexp(line,[name '\(([-+]*\w+)\)'],'tokens','once');
0283             % set argument as expression (see dsWriteDynaSimSolver() for usage as such)
0284             if ~isempty(arg)
0285               rhs=arg;
0286             end
0287           end
0288           
0289           % convert into valid monitor name
0290           name=strrep(name,'.','_'); % index sub-namespace (monitor Na.I)
0291           
0292           if ~isempty(rhs), expression=rhs{1}; else expression=[]; end
0293           
0294           model.monitors(1).([namespace name]) = expression;
0295           name_map(end+1,:) = {name,[namespace name],namespace,'monitors'};
0296           
0297           if ~isempty(comment)
0298             model.comments{end+1}=sprintf('%s (monitor): %s',[namespace name],comment);
0299           end
0300         end
0301       
0302       end
0303       
0304     case 'conditional'      % if(conditions)(actions)
0305       groups=regexp(line,'\)\(','split');
0306       condition=regexp(groups{1},'^if\s*\((.*)','tokens','once');
0307       if length(groups)==2
0308         if groups{2}(end)==')'
0309           groups{2}=groups{2}(1:end-1);
0310         end
0311         
0312         then_action=groups{2};
0313         else_action=[];
0314       elseif numel(groups==3)
0315         if groups{3}(end)==')'
0316           groups{3}=groups{3}(1:end-1);
0317         end
0318         
0319         then_action=groups{2};
0320         else_action=groups{3};
0321       end
0322       
0323       model.conditionals(end+1).namespace=namespace;
0324       model.conditionals(end).condition=condition{1};
0325       model.conditionals(end).action=strrep(then_action,',',';'); % restore semicolon-delimited multiple actions like if(x>1)(x=0;y=0)
0326       
0327       if length(groups)>2
0328         model.conditionals(end).else=else_action;
0329       else
0330         model.conditionals(end).else=[];
0331       end
0332       
0333       if ~isempty(comment)
0334         model.comments{end+1}=sprintf('%s conditional(%s): %s',namespace,condition,comment);
0335       end
0336       
0337     case 'linker'           % [link ]? target operation expression (e.g., link target += f(x))
0338       % viable options: ((\+=)|(-=)|(\*=)|(/=)|(=>))
0339       line=regexprep(line,'^link (\s*\w)','$1');
0340       lhs=regexp(line,'^([^\+\-*/=]+)','tokens','once'); % +=
0341       rhs=regexp(line,'=>?(.+)$','tokens','once');
0342       model.linkers(end+1).namespace=namespace;
0343       
0344       if ~isempty(lhs), target=strtrim(lhs{1}); else target=[]; end
0345       
0346       if ~isempty(rhs), expression=strtrim(rhs{1}); else expression=[]; end
0347       
0348       if expression(end)==';', expression=expression(1:end-1); end
0349       
0350       if isempty(target), target=expression; end % for sharing state var across mechanisms in same population
0351       
0352       if isempty(expression), expression=target; end
0353       
0354       model.linkers(end).target=target;
0355       model.linkers(end).expression=expression;
0356       model.linkers(end).operation='+=';
0357       
0358       if ~isempty(comment)
0359         model.comments{end+1}=sprintf('%s linkers(%s->%s): %s',namespace,target,expression,comment);
0360       end
0361       
0362       if ~isempty(comment)
0363         model.linkers(end).comment=comment;
0364       end
0365     otherwise
0366       warning('ignoring line, failed to classify :%s',line);
0367   end
0368 end
0369 
0370 %% auto_gen_test_data_flag argout
0371 if options.auto_gen_test_data_flag
0372   argout = {model, name_map}; % specific to this function
0373   
0374   dsUnitSaveAutoGenTestData(argin, argout);
0375 end
0376 
0377 end % main fn
0378 
0379 
0380 %% local functions
0381 function [line,comment]=remove_comment(line)
0382 % purpose: split line into model content and comment
0383 index=find(line=='%',1,'first');
0384 if isempty(index)
0385   % check other valid comment delimiters
0386   index=find(line=='#',1,'first');
0387 end
0388 
0389 if isempty(index) % no comment found
0390   comment='';
0391 else % split line into comment and non-comment line sections
0392   comment=line(index:end);
0393   line=line(1:index-1);
0394 end
0395 
0396 end

Generated on Tue 12-Dec-2017 11:32:10 by m2html © 2005