boussinesq_lookup.m

function sigma_z=boussinesq_lookup(q,a,b,z)
  % function that determines stress under corner of an a by b rectangular platform
  % z-meters below the platform. The calculated solutions are in the fmn data
  % m=fmn(:,1)
  % in column 2, fmn(:,2), n=1.2
  % in column 3, fmn(:,2), n=1.4
  % in column 4, fmn(:,2), n=1.6

  fmn= [0.1,0.02926,0.03007,0.03058
        0.2,0.05733,0.05894,0.05994
        0.3,0.08323,0.08561,0.08709
        0.4,0.10631,0.10941,0.11135
        0.5,0.12626,0.13003,0.13241
        0.6,0.14309,0.14749,0.15027
        0.7,0.15703,0.16199,0.16515
        0.8,0.16843,0.17389,0.17739];

  m=a/z;
  n=b/z;
  if n < 1.3
      f=fmn(:,2);
    elseif n > 1.5
        f=fmn(:,4);
    else
        f=fmn(:,3);
  end
  [~,i_fit]=sort(abs(m-fmn(:,1)));
  x=fmn(i_fit(1:4),1);
  y=f(i_fit(1:4));
  f_out = Newtint(x,y,m);
  sigma_z=q*f_out;
end