function [sig_v, sig_h, sig_H] = Sibson_NormalFaultTensor(d,lm, mu_s)

%determination of stress tensor based on R.H. Sibson, Frictional constraints on thrust, wrench, and normal faults, Nature (1974)
%note: it is assumed that sig_H= mean normal stress = intermediate principal stress = (sig_h+sig_v)/2
%theta = (angle between sigma 1 and fault)
%mu_s = (static coefficient of friction)
theta_p = 0.5*atan2(1, mu_s);%optimally oriented fault angle
rho = 2650; %mass density of rocks (kg/m3)
%lm = (pore fluid factor - ratio of hydrostatic to lithostatic pressure gradient)
g = 10; %acceleration due to gravity (m/s)

R = (sin(2*theta_p)+mu_s*(cos(2*theta_p)+1))/(sin(2*theta_p)+mu_s*(cos(2*theta_p)-1));
sig_v = rho*g*d*(1-lm)*10^-6; %vertical princ stress (MPa) - *compression positive
sig_h = sig_v-(R/(1+R))*sig_v; %least compressive horizontal princ stress (MPa) - *compression positive
sig_H = (sig_v+sig_h)/2; %greatest compressive horizontal princ stress (MPa) - *compression positive