MCS_square_footing.m

%% Monte Carlo method for estimating the prob of failure of a square footing
%{
---------------------------------------------------------------------------
*Created by:                        Date:             Comment:
Felipe Uribe-Castillo               November 2009     Final work
*Mails:                             Manizales
[e-mail]             
[e-mail]   
*University:
National University of Colombia     
---------------------------------------------------------------------------
Failure probability analysis of a superficial square footing on plastic 
soil of stiff consistency, using crude Monte Carlo Simulation.
---------------------------------------------------------------------------
Based on:
1. C. Venkatramaiah. "Geotechnical engineering". 1993. John Wiley and Sons.
2. R. Rackwitz. "Reviewing probabilistic soils modelling". In Computers and 
   Geotechnics 26(3-4), 2000, pp 199 - 223.
---------------------------------------------------------------------------
%}
clear; close all; format long g; clc;

%% Initial Data
N    = 1e6;        % Number of Random Numbers
H    = 3;          % Depth of strata, m
Df   = 2;          % "Profundidad de desplante", m (D<=B)
B    = 2.50;       % Width of the footing, m
Qapp = 1000;       % Applied load, kN (Supossed)
Smax = 0.30;       % Max seetlement, m (NSR)
z    = 3;          % Depth in which is measured the settlement, m
mv   = 0.000125;   % Coefficient of volume compressibility , m^2/kN (Supossed)

% Laboratory results for friction angle. 
% (Mean values taken of Rackwitz's paper, pag17)
m1     = 17.45;         % deg
v1     = 4.57^2;        % deg
mu1    = log((m1^2)/sqrt(v1+m1^2));
sigma1 = sqrt(log(v1/(m1^2)+1));

% Laboratory results for cohesion in consolidated soil. 
% (Mean values taken of Rackwitz's paper, pag17)
mu2    = 10;            % kN/m^2
sigma2 = 4.5;           % kN/m^2

% Laboratory results for specifics weigths. 
% (Mean values taken of Rackwitz's paper, pag17)
mu3    = 18;            % kN/m^3
sigma3 = 1;             % kN/m^3

%% Calculate the soil parameters
phi    = lognrnd(mu1,sigma1,N,1)*pi/180;     % rad
c      = normrnd(mu2,sigma2,N,1);            % kN/m^2
gamma  = normrnd(mu3,sigma3,N,1);            % kN/m^3

%% Bearing capacity factors asignation
Nq = zeros(N,1);
Nc = zeros(N,1);
Ng = zeros(N,1);
% Using Terzaghi's theory
for i = 1:N
    Nq(i) = (exp( ((3*pi/4)-(phi(i)/2))*tan(phi(i)) ))^2/...
            ( 2*(cos((pi/4)+(phi(i)/2)))^2 );
    Nc(i) = cot(phi(i))*(Nq(i)-1);
    Ng(i) = 0.5*tan(phi(i)) * ( (((tan((pi/4)+(phi(i)/2)))^2)/...
            ((cos(phi(i)))^2)) - 1 );
end

%% Calculate the bearing capacity of the foundation using Terzaghi's Method
qc = zeros(N,1);
% Terzaghi's equation for square footing
for i = 1:N
    qc(i) = 1.3*c(i)*Nc(i) + gamma(i)*Df*Nq(i) + 0.4*B*gamma(i)*Ng(i);  % kN/m^2
end
% Calculate the reliability analysis function for the applied loads
g    = qc-(Qapp/(B^2));
% Failure probability respect to the loads
pf_l     = mean(g<=0);
var_pf_l = pf_l*(1-pf_l)/N;
std_pf_l = sqrt(var_pf_l);
fprintf('\n Probability of failure of the square footing by load: %8.6f \n',...
         pf_l)

%% Calculate the settlement of the foundation using Terzaghi's Method
sigmaz = zeros(N,1);
deltaH = zeros(N,1);
% Using Boussinesq equation
m = B/z;
I = ((2*m^2*sqrt(2*m^2+1))/(m^4+2*m^2+1))*((2*m^2+2)/(2*m^2+1))*...
    atan2((-m^4+2*m^2+1),(2*m^2*sqrt(2*m^2+1)));
for i = 1:N
    sigmaz(i) = (qc(i)/4*pi)*I;
    deltaH(i) = mv*sigmaz(i)*H;   % m
end
% Calculate the reliability analysis function for the settlement
g2 = Smax-deltaH;
% Failure probability respect to the settlement
pf_s     = mean(g2<=0);
var_pf_s = pf_s*(1-pf_s)/N;
std_pf_s = sqrt(var_pf_s);
fprintf('\n Probability of failure of the square footing by settlement: %8.8f \n\n',...
        pf_s)

%End