Description of genericDerivative

0001 function derivativeStruct = genericDerivative(solverState, time, ...
0002     pbeDefinition)
0003 % compute derivative of states from structures
0004 %
0005 % THIS IS NO USER FUCNTION
0006 
0007 % The elk-library: convex geometry applied to crystallization modeling.
0008 %   Copyright (C) 2012 Alexander Reinhold
0009 %
0010 % This program is free software: you can redistribute it and/or modify it
0011 %   under the terms of the GNU General Public License as published by the
0012 %   Free Software Foundation, either version 3 of the License, or (at your
0013 %   option) any later version.
0014 %
0015 % This program is distributed in the hope that it will be useful, but
0016 %   WITHOUT ANY WARRANTY; without even the implied warranty of
0017 %   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
0018 %   General Public License for more details.
0019 %
0020 % You should have received a copy of the GNU General Public License along
0021 %   with this program.  If not, see <http://www.gnu.org/licenses/>.
0022 
0023 %% verbosity
0024 % if verbose
0025 %     time/pbeDefinition.timeSpan(end)*100
0026 % end
0027 
0028 %! DEV: dbugging pbe lactose case
0029 % if time > 6000
0030 %     warning('time is up')
0031 % end
0032 %/
0033 
0034 %% initialize
0035 time = time
0036 nDecCoordinate = pbeDefinition.nDecCoordinate;
0037 nExtraCoordinate = pbeDefinition.nExtraCoordinate;
0038 
0039 %% ensure hcVectors in cone and proper upVector
0040 % Note that pivotUnifiedPartitionVector can contain -1 when solver steps
0041 %   have become too large. Fore measure calculation, this is corrected to
0042 %   arbitrary partition=1, the growth rate will then be set to 0, a warning
0043 %   displayed and the soler should correct the stepsize cause of the sudden
0044 %   change in growth.
0045 [pivotPositionMeasureMatrix pivotUnifiedPartitionVector fractionCorrected] = ...
0046     ensureProperDecData(solverState.hcMatrix, solverState.upVector, ...
0047       pbeDefinition);
0048 
0049 %% pivot properties
0050 pivotPropertyMatrix = pbeDefinition.pivotProperty(pivotPositionMeasureMatrix, ...
0051     max(pivotUnifiedPartitionVector, 1));
0052 
0053 %% bulk properties
0054 bulkPropertyVector = pbeDefinition.bulkProperty(pivotPropertyMatrix, ...
0055     solverState.numberVector, solverState.bulkStateVector, time);
0056 
0057 %% physical rates
0058 bulkRate = pbeDefinition.bulkRate(bulkPropertyVector, time);
0059 nucleationRate = pbeDefinition.nucleationRate(bulkPropertyVector, time);
0060 growthRateMatrix = pbeDefinition.growthRate(pivotPositionMeasureMatrix, ...
0061     pivotPropertyMatrix, bulkPropertyVector, time);
0062 
0063 %% adjust rates to boundary
0064 % The growthRateMatrix might be given in the proper embedding for improper
0065 %   representations (users choice). This detail is handled by
0066 %   limitGrowthRate to avoid uneccesary projections during this algorithm.
0067 if pbeDefinition.optionStruct.limitGrowthRate
0068     positionDerivativeMatrix = [...
0069             limitGrowthRate(growthRateMatrix(1:(end-nExtraCoordinate), :), ...
0070             solverState.hcMatrix(1:nDecCoordinate, :), ...
0071             pivotUnifiedPartitionVector, ...
0072             pbeDefinition.dec, ...
0073             pbeDefinition.optionStruct.distanceLimit);...
0074         growthRateMatrix((end-nExtraCoordinate+1):end, :)];
0075 else
0076     positionDerivativeMatrix = growthRateMatrix;
0077 end
0078 
0079 %% assign derivatives as solver struct
0080 derivativeStruct.bulkStateVector = bulkRate;
0081 derivativeStruct.hcMatrix = positionDerivativeMatrix;
0082 derivativeStruct.flagDensity = solverState.flagDensity;
0083 if ~isempty(positionDerivativeMatrix)
0084     derivativeStruct.numberVector = positionDerivativeMatrix(1,:)*0;
0085 else
0086     derivativeStruct.numberVector = zeros(1, 0);
0087 end
0088 if solverState.flagNucleation
0089     derivativeStruct.numberVector(1) = nucleationRate;
0090     derivativeStruct.hcMatrix(:,1) = 0.5*derivativeStruct.hcMatrix(:,1);
0091 end
0092 
0093 %% rate of probability density
0094 % to provide the rate of the probability density, we need to calculate the
0095 %   divergence of the growth rate.
0096 if solverState.flagDensity
0097 %     disp('DENSITY')
0098     divergenceVector = zeros(1, length(solverState.probabilityVector));
0099     for iDim = 1:size(solverState.hcMatrix, 1)
0100         % get a state somewhat distant in h-space
0101         forwardSolverState = solverState;
0102         deltaVector = min(...
0103             pbeDefinition.optionStruct.gradientRelativeDelta * ...
0104               forwardSolverState.hcMatrix(iDim, :), ...
0105             pbeDefinition.optionStruct.gradientAbsoluteDelta + ...
0106               0*forwardSolverState.hcMatrix(iDim, :));
0107                
0108         forwardSolverState.hcMatrix(iDim, :) = ...
0109             forwardSolverState.hcMatrix(iDim, :) + deltaVector;
0110         % remove density flag for that we only need the growth rate
0111         forwardSolverState.flagDensity = 0;
0112         forwardDerivativeStruct = genericDerivative(forwardSolverState, time, ...
0113             pbeDefinition);
0114         
0115         % accumulate divergence
0116         divergenceVector = divergenceVector + ...
0117             (forwardDerivativeStruct.hcMatrix(iDim, :) - ...
0118              derivativeStruct.hcMatrix(iDim, :)...
0119              ) ./ deltaVector;
0120     end
0121     derivativeStruct.probabilityVector = -1*...
0122         solverState.probabilityVector .* divergenceVector;
0123 end
0124 
0125 %% asign info for unified partition updating
0126 derivativeStruct.fractionCorrected = fractionCorrected;
0127 derivativeStruct.partitionIndexVector = pivotUnifiedPartitionVector;
0128 
0129 end
genericDerivative

PURPOSE

SYNOPSIS

DESCRIPTION

CROSS-REFERENCE INFORMATION

SOURCE CODE
