function pbeDefinition = growthDensityCommon

%% Name and information
pbeDefinition.info = ['seeded growth (one initial population), ' ...
    'using a growth field that accelerates, decellerates ' ...
    'and goes against a boundary'];
pbeDefinition.source = 'generic';

%% Crystal setup
% crystal
pbeDefinition.cdef = obtainCrystal('generic_potash_alum', 1);

% decomposition data
pbeDefinition.dec = obtainDec(pbeDefinition.cdef.A, 'mapping', ...
    pbeDefinition.cdef.groupMappingMatrix, 'verbose', 0);
pbeDefinition.dec = addMeasureDataDec(pbeDefinition.dec, 'volume');
pbeDefinition.dec = addMeasureDataDec(pbeDefinition.dec, 'surface');
pbeDefinition.dec = addMeasureDataDec(pbeDefinition.dec, 'meanWidth');

%% Initial distribution
pbeDefinition.seeding.time = 0; % seconds
pbeDefinition.seeding.mean1 = [10 10]'*1e-7;
pbeDefinition.seeding.std1 = 0.5*[2 2]'*1e-7;   % m
pbeDefinition.seeding.totalValue1 = 1e10;
pbeDefinition.seeding.mean2 = [12.5 12.5]'*1e-7;
pbeDefinition.seeding.std2 = 0.5*[3 3]'*1e-7;   % m
pbeDefinition.seeding.totalValue2 = 1e10;
pbeDefinition.seeding.nSample = 100;
% the particle distribution
particleDisCOne = obtainDistributionContinuous('normal', ...
    'mean', pbeDefinition.seeding.mean1, ...  %m
    'cov', diag(pbeDefinition.seeding.std1.^2), ...
    'totalValue', pbeDefinition.seeding.totalValue1);
particleDisCTwo = obtainDistributionContinuous('normal', ...
    'mean', pbeDefinition.seeding.mean2, ...  %m
    'cov', diag(pbeDefinition.seeding.std2.^2), ...
    'totalValue', pbeDefinition.seeding.totalValue2);
% what I need here is a function that sets the value artificially to zero,
%   if the position is invalid. How to link that in?
%     * logicalOut = inCone(cone.A, hCMatrix) - easily applied for above..
%         can be applied in <distribution>, <pbeSolver>, <dec>. <dec> it
%         should be.
pbeDefinition.seeding.fhandle = @(h)(...
    isValidDec(h, pbeDefinition.dec, 0) .* (...
    particleDisCOne.fhandle(h) + ...
    particleDisCTwo.fhandle(h) ));

%% growth and properties
% shortcut
dec = pbeDefinition.dec;

% growth rate
% crystals start at mean [1 1]*1e-5 and grow with [1 1]*5e-9.
% For hC1 > 2e-5 the rate in hC1 direction doubles until 3e-5.
% For hC2 >
% 2.5e-9, the rate halfs down until 4e-5.
lb1 = 4e-5;ub1 = 5e-5; r1 = ub1 - lb1;
factorOne = @(H)(1.0 + ...
    (H(1,:) > lb1).*(H(1,:) < ub1).*((H(1,:)-lb1)/r1)*2.0 + ...
    (H(1,:) >= ub1)*2.0...
    );
% factorOne(pbeDefinition.seeding.positionMatrix)
lb2 = 2.5e-5;ub2 = 3.5e-5; r2 = ub2 - lb2;
factorTwo = @(H)(1.0 - ...
    (H(2,:) > lb2).*(H(2,:) < ub2).*((H(2,:)-lb2)/r2)*0.5 - ...
    (H(2,:) >= ub2)*0.5...
    );
% factorOne(pbeDefinition.seeding.positionMatrix)

pbeDefinition.growthRate = @(H, P, y, t)(...
    [6; 10]*1e-9 * ones(1, size(H, 2)) .* ...
    [factorOne(H); ones(1, size(H, 2))] .* ...
    [ones(1, size(H, 2)); factorTwo(H)] );

% properties
pbeDefinition.pivotProperty = {'volume'};
pbeDefinition.bulkProperty = @(P, N, x, t)(P*N');

% time span for solution
% pbeDefinition.timeSpan = 0:30:3.5e3;          %s
% pbeDefinition.timeSpan = 0:60:2e4;          %s
pbeDefinition.timeSpan = 0:30:1.3e4;          %s

% options
pbeDefinition.optionStruct.distanceLimit = 1e-7; %m
pbeDefinition.optionStruct.density = 1;
pbeDefinition.optionStruct.relTol = 1e-5;
pbeDefinition.optionStruct.absTolPosition = 1e-10;
pbeDefinition.optionStruct.absTolDensity = 1e8;
pbeDefinition.optionStruct.odeSolver = @ode45;
pbeDefinition.optionStruct.gradientRelativeDelta = 1e-3;
pbeDefinition.optionStruct.gradientAbsoluteDelta = 1e-10;