Skip to content

Commit

Permalink
small changes in example Groningen
Browse files Browse the repository at this point in the history
  • Loading branch information
Buijze committed Sep 5, 2024
2 parents 5d84564 + 2f9c01e commit 9edc7b0
Show file tree
Hide file tree
Showing 6 changed files with 60 additions and 26 deletions.
58 changes: 41 additions & 17 deletions bin/MultiFaultCalculator.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -124,38 +124,62 @@
end

function self = set_run_setting(self, setting_name, setting_value)
% specify run settings per pillar. (exc load table, y, input)
% Specify run settings per pillar.
% INPUT
% setting name
% setting value cell array, or array of floats or single
% cell/single float
% setting_name - Name of the setting to be applied
% setting_value - Cell array, array of floats, single cell, single float, or string

[is_valid, value_type] = self.is_valid_setting_name(setting_name);

if is_valid
valid_value = 1; % Initialize as valid

for i = 1 : length(self.pillars)
valid_value = 1;
if iscell(setting_value) && (length(setting_value) == length(self.pillars))
assign_value = setting_value{i};
elseif iscell(setting_value) && (length(setting_value) == 1)
assign_value = setting_value{1};
elseif isfloat(setting_value) && (length(setting_value) == 1)
assign_value = setting_value(1);
elseif isfloat(setting_value) && (length(setting_value) == self.n_pillars)
assign_value = setting_value(i);
if iscell(setting_value)
% Case 1: Cell array with the same length as pillars
if length(setting_value) == length(self.pillars)
assign_value = setting_value{i};
% Case 2: Single cell element (could be string, char array, or numeric)
elseif length(setting_value) == 1
assign_value = setting_value{1};
else
valid_value = 0; % Invalid case: cell array with incorrect length
end

elseif isnumeric(setting_value)
% Case 3: Single numeric value
if length(setting_value) == 1
assign_value = setting_value(1);
% Case 4: Numeric array with the same length as pillars
elseif length(setting_value) == self.n_pillars
assign_value = setting_value(i);
else
valid_value = 0; % Invalid case: numeric array with incorrect length
end

elseif ischar(setting_value) || isstring(setting_value)
% Case 5: String or character array
assign_value = setting_value;

else
valid_value = 0;
valid_value = 0; % Invalid case: unsupported type
end
if valid_value & strcmp(value_type,class(assign_value))

% Assign value if valid and type matches the expected type
if valid_value && strcmp(value_type, class(assign_value))
self.pillars{i}.(setting_name) = assign_value;
end
end
if ~valid_value || ~strcmp(value_type,class(assign_value))

% Display a message if the value was not valid or the type didn't match
if ~valid_value || ~strcmp(value_type, class(assign_value))
disp('Specified setting type does not seem the right type or dimension, check');
end
end

end



function [nuc_load_step] = get_minimum_nucleation_load_step(self)
if self.run_done
nuc_load_step = min(self.result_summary.nucleation_load_step);
Expand Down
3 changes: 1 addition & 2 deletions bin/panther.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,7 @@
% of seconds. If it is already running (type gcp to check), parfor will initiate much faster
parfor (i = 1 : n_members, matlab_workers)
% disp([num2str(i),'/', num2str(n_members)]);
L{i} = y./sin(ensemble{i}.dip*pi/180);
L{i} = y./sin(ensemble{i}.dip*pi/180);

% initial stress
initial_stress{i} = InitialStress(y, ensemble{i});
Expand Down Expand Up @@ -128,7 +128,6 @@

reduced_load_table = analysis.load_table(indices_for_saving,:);



for i = 1 : n_members
% save stresses, pressure, slip
Expand Down
4 changes: 2 additions & 2 deletions examples/example_depth_dependent_property.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -77,8 +77,8 @@
% set depth varying dip
sc2_variable_dip.input_parameters.dip.uniform_with_depth = 0; % make dip depth-variable
y = sc2_variable_dip.y;
sc2_variable_dip.input_parameters.dip.value_with_depth = 90*ones(size(y));
i_mid = floor(length(y)/2);
sc2_variable_dip.input_parameters.dip.value_with_depth = 80*ones(size(y));
i_mid = floor(length(y)/2)
sc2_variable_dip.input_parameters.dip.value_with_depth(i_mid:end) = 60;
% sc2_variable_dip.input_parameters.dip.value_with_depth = 60 + y*dip_gradient_per_m;

Expand Down
16 changes: 12 additions & 4 deletions src/calculator/FaultSlip.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -104,7 +104,7 @@
% length(time).
slip_zone_indices = self.get_slip_zone_indices(slipping);
slip_zone_length = nan(size(slip_zone_indices)); % along-fault length of slip zones
nucleation_length = nan(size(slip_zone_indices)); % nucleation length, per slip zone
nucleation_length = nan(size(slip_zone_indices)); % nucleation length, per slip zone
% iterate over time steps
for i = 1 : size(slip_zone_indices,2)
% iterate over number of slip zones
Expand All @@ -121,6 +121,7 @@
slip_zone_length(j,i) = L_start - L_end;
sne_slip = (sne(slip_zone_indices{j,i}, i));
tau_slip = (tau(slip_zone_indices{j,i}, i));

if length(f_s) == size(tau,1)
% add ,i here if friction changes per timestep
f_s_slip = f_s(slip_zone_indices{j,i}) ;
Expand All @@ -137,9 +138,15 @@
else
d_c_slip = d_c;
end
if length(mu_II) == size(tau,1)
mu_II_slip = mean(mu_II(slip_zone_indices{j,i}));
else
mu_II_slip = mu_II;
end
delta_tau = mean(sne_slip .* (f_s_slip - f_d_slip));
tau_0_d = mean(tau_slip - sne_slip .* f_d_slip);
nucleation_length(j,i) = self.calculate_nucleation_length( delta_tau, tau_0_d, d_c_slip, mu_II, nuc_crit, nuc_len_fixed);

nucleation_length(j,i) = self.calculate_nucleation_length(delta_tau, tau_0_d, d_c_slip, mu_II_slip, nuc_crit, nuc_len_fixed);
%nucleation_length(j,i) = 1.158 * mu_II * d_c./((f_s - f_d) * average_sne_in_slip_zone);
end
end
Expand Down Expand Up @@ -217,8 +224,9 @@
% mu_II mode II shear stiffness
nx = length(L);
Kline = -nx/(2*pi*(L(1)-L(end))) ./ ( [0:nx-1]'.^2-0.25 );
Ko = toeplitz(Kline);
K = mu_II*Ko;
Ko = toeplitz(Kline);
% K = mu_II.*Ko;
K = Ko*mu_II;
K(and(K<0,K>(min(min(K)/10000)))) = 0; % set very small changes to 0 to avoid continued interactions of the two peaks to
end

Expand Down
3 changes: 3 additions & 0 deletions src/calculator/FaultStress.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,9 @@
self.tau = initial.tau0 + change.dtau;
end

function self = get_reactivation_nucleation_stress(self, reactivation_load_step, nucleation_load_step)
[self.sne_reac, self.tau_reac] = self.get_stress_at_load_step(reactivation_load_step);
end

function self = get_nucleation_stress(self, nucleation_load_step)
% INPUT
Expand Down
2 changes: 1 addition & 1 deletion src/calculator/FaultStressChange.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -123,7 +123,7 @@
end
if length(gamma_PT) == 1
gamma = gamma_PT;
elseif length(gamma_P) == length(y)
elseif length(gamma_PT) == length(y)
gamma = gamma_PT(j);
end
[dsn_j, dtau_j] = self.get_stress_change_component(GF{j}, dPT_FW, dPT_HW, gamma);
Expand Down

0 comments on commit 9edc7b0

Please sign in to comment.