www.gusucode.com > simbio 案例源码程序 matlab代码 > simbio/discontSimBiologyModel.m
%% Deterministic Simulation of a Model Containing a Discontinuity % % This example shows how to correctly build a SimBiology(R) model % that contains discontinuities. % Copyright 2010-2014 The MathWorks, Inc. %% Background % The model you create in this example simulates the first-order elimination % of a protein that is produced at a specified rate. The production rate % contains two discontinuities. To simulate the model accurately, you must % create events that are triggered at the time of the discontinuity. % % The production rate has three "modes" of production, as illustrated in % the following plot: % plot([0 3 3 6 6 10], [5 5 3 3 0 0]); ylim([-0.5 5.5]); xlabel('Time'); ylabel('Rate'); title('Discontinuous Protein Production Rate'); %% % Initially ("Mode 1"), the production rate is a constant value of 5. From % 3 to 6 seconds ("Mode 2"), the production rate is 3. After 6 seconds % ("Mode 3"), the production rate is 0. These production rates are % implemented in a MATLAB function discontSimBiologyRateFunction.m, which % requires two arguments, simulation time and production mode. % % In this example, you will add events to the model to change the mode of % protein production. This approach ensures that discontinuities in the % model occur only via events, which further ensures that the ODE solver % accurately calculates the numerical behavior near the discontinuities. % % Note that to simulate a model accurately you must use events to handle % _any_ discontinuity, whether related to function values or their % derivatives. %% Construct the Model, Compartment, and Species model = sbiomodel('discontinuous rate'); central = addcompartment(model,'Central'); addspecies(central,'protein') %% Construct the Reaction for First-Order Elimination reaction1 = addreaction(model,'protein -> null') ke = addparameter(model,'ke', 0.5); kineticLaw1 = addkineticlaw(reaction1,'MassAction'); kineticLaw1.ParameterVariableNames = {ke.Name}; reaction1.ReactionRate; %% Construct the Events That Are Triggered at the Time of Discontinuities % These events set a parameter mode that controls the mode of protein % production. The mode is initially 1, changes to 2 at time 3, and changes % to 3 at time 6. counter = addparameter(model,'mode', 1, 'ConstantValue', false); addevent(model,'time > 3', 'mode = 2') addevent(model,'time > 6', 'mode = 3') %% Construct the Reaction for Protein Production % We assign this rate to a parameter using a repeated assignment rule. This % lets us store the production rate in the simulation results. reaction2 = addreaction(model, 'null -> protein'); rate2 = addparameter(model,'rate2', 0, 'ConstantValue', false); reaction2.ReactionRate = 'rate2' addrule(model,'rate2 = discontSimBiologyRateFunction(time, mode)', 'repeatedAssignment') %% View the Contents of discontSimBiologyRateFunction type discontSimBiologyRateFunction %% Simulate and Plot the Model model sd = sbiosimulate(model); plot(sd.Time, sd.Data); ylim([-0.5 8]); xlabel('Time'); ylabel('State'); title('Simulation Results'); legend(sd.DataNames); %% Conclusion % This example illustrates how to create a SimBiology model that contains % discontinuities. It illustrates how to add events to the model to address % the discontinuities, so you can simulate the model accurately.