// JavaScript Document

// layer consructor
function layer(thick, CPTqc, N60, gammaDry, gammaWet, E) {
	this.thickness = Number(thick);
	this.qc = Number(CPTqc);
	this.N = Number(N60);
	this.unitWeightDry = Number(gammaDry);
	this.unitWeightWet = Number(gammaWet);
	this.Youngs = Number(E);
}

// function to subdivide layers based on input B.  Want 10 layers per B
// return a new array of layers
function subdivideLayers(layers, B)
{
	subThickness = 0.1*B;
	newLayers = [];
	for (i=0; i<layers.length; i++) {
		subLayers = parseInt(layers[i].thickness / subThickness);
		if (subLayers < 1)  subLayers = 1;
		_thickNew = layers[i].thickness / subLayers;
		_qc = layers[i].qc;
		_E = layers[i].Youngs;
		_N = layers[i].N;
		_gammaDry = layers[i].unitWeightDry;
		_gammaWet = layers[i].unitWeightWet;
		
		for (j=0; j<subLayers; j++) {
		  newLayers.push(new layer(_thickNew, _qc, _N, _gammaDry, _gammaWet, _E));
		}
	}
	return newLayers;
	
}

// Schmert - convert qc to Es
function convertToEs(qc, footing, useModified) {
  if (useModified) {
	 	// axisymmetric if width = length, Es = 2.5*qc
		// plane strain if length = 10 X width, Es = 3.5*qc
		ratio = footing.length/footing.width;
		if (ratio < 1.0) return (2.5*qc);
		else if (ratio > 10.0) return (3.5*qc);
		
		return ((ratio-1.0)/9.0*3.5*qc + (10.0-ratio)/9.0*2.5*qc);
  }
  else  return (2.0*qc);
}

// input is depth from surface
function getOverburdenStress(depth, piezo, layers, gammaWater) {
  	_top = 0.0;
	_totalStress = 0.0;
	for (i=0; i<layers.length; i++) {
		// + by default concatenates strings.  Must convert to numbers
		_bottom = _top + layers[i].thickness;
		
		// hack - if it is the last layer, and we haven't reached depth, set bottom to depth
		if (i == layers.length-1) {
			if (_bottom < depth) _bottom = depth;
		}
		if (_bottom > depth) _bottom = depth;
		
		dryThickness = _bottom - _top;
		if (piezo < _top) dryThickness = 0.0;
		else if (piezo < _bottom) dryThickness = piezo - _top;
		wetThickness = (_bottom - _top) - dryThickness;
		_totalStress += dryThickness*layers[i].unitWeightDry + wetThickness*layers[i].unitWeightWet;
		_top = _bottom;
		if (_top >= depth) break;
	}
	if (piezo < depth) _totalStress -= gammaWater*(depth - piezo);
	return _totalStress;
}

// Schmert - this returns depth below bottom of load - not depth from surface
function getIzPeakDepth(footing, useModified) {
	if (useModified) {
		// axisymmetric if width = length, depth = B/2
		// plane strain if length = 10 X width, depth = B
		ratio = footing.length/footing.width;
		if (ratio < 1.0) return (0.5*footing.width);
		else if (ratio > 10.0) return (footing.width);
		
		return ((ratio-1.0)/9.0*footing.width + (10.0-ratio)/9.0*0.5*footing.width);
	}
	else return (0.5*footing.width);
	
}

function getIzBottomDepth(footing, useModified) {
  _bottom = 2.0*footing.width;
  if (useModified) {	
  	 ratio = footing.length/footing.width;
	 if (ratio < 1.0) ratio = 1.0;
	 else if (ratio > 10.0) ratio = 10.0;
	 
	 _bottom = (ratio-1.0)/9.0*4.0*footing.width + (10.0-ratio)/9.0*2.0*footing.width;
  }

  return _bottom;
	
}

// for Schmert
function getIz(depthBelowFooting, footing, IzPeak, IzDepth, useModified) {
  if (useModified) {	
  	 ratio = footing.length/footing.width;
	 if (ratio < 1.0) ratio = 1.0;
	 else if (ratio > 10.0) ratio = 10.0;
	 
	 _bottom = (ratio-1.0)/9.0*4.0*footing.width + (10.0-ratio)/9.0*2.0*footing.width;
	 
     if (depthBelowFooting < IzDepth) {
		 _b = (ratio-1.0)/9.0*0.2 + (10.0-ratio)/9.0*0.1;
		 _m = (IzPeak - _b)/IzDepth;
		 return (_m*depthBelowFooting+_b);
	 }
	 else if (depthBelowFooting < _bottom) {
		 _m = -IzPeak/(_bottom - IzDepth);
		 _b = -_m*_bottom;
		 return (_m*depthBelowFooting+_b);
	 }
	 else return 0.0;
  }
  else {
	  if (depthBelowFooting < 0.5*footing.width) return (1.2*depthBelowFooting/footing.width);
	  else if (depthBelowFooting < 2.0*footing.width) return (-0.4/footing.width*depthBelowFooting+0.8);
	  else return 0.0; 
  }
	
	
}

// INPUT: footing object
//			piezo depth from surface
//			time
//			an array of layer objects
//			boolean to indicate if using modified (1978)
//			unit weight of water
//			boolean to indicate if layers should be subdivided
// OUTPUT:  a 2D array of depths and settlements
function computeSchmertSettlement(footing, piezo, time, layers, useModified, gammaWater, subdivide)
{
	// first align layers at the bottom
	schmertBottom = getIzBottomDepth(footing,useModified)  + Number(footing.depth);
	_bottom = 0;
	for (i=0; i<layers.length; i++) {
		_bottom += layers[i].thickness;
		if (_bottom > schmertBottom) {
			layers[i].thickness -= _bottom - schmertBottom;
			break;
		}
		
	}
	
	// convert qc to Es
	for (i=0; i<layers.length; i++) {
		if (layers[i].Youngs == 0)  layers[i].Youngs = convertToEs(layers[i].qc,footing,useModified);
	}
	
	sz = getOverburdenStress(footing.depth,piezo,layers,gammaWater); // overburden pressure at base of footing
	delta_q = footing.magnitude - sz; // new foundation pressure
	
	if (delta_q < 0) {
		alert('Warning: Load is less than overburden stress at bottom of footing');
		return 0;
	}
	C1 = 1.0 - 0.5*sz/delta_q;
	if (C1 < 0.5) C1 = 0.5;
	C2 = 1.0; 
	if (time > 0) C2 += 0.2*Math.log(10*time)/Math.log(10);
	
	IzPeak = 0.6;
	IzDepth = getIzPeakDepth(footing,useModified);
	if (useModified) {
		sz = getOverburdenStress(IzDepth+footing.depth,piezo,layers, gammaWater);
		IzPeak = 0.5+0.1*Math.sqrt(delta_q/sz);
	}
	
	// now subdivide layers
	if (subdivide) {
	  dividedLayers = subdivideLayers(layers,footing.width);
	}
	else {
		dividedLayers = layers;
	}
	
	// make 2D array of compaction in each layer
	// convert this to total settlement below
	compactionArray = [];
	cc = C1*C2*delta_q;
	_top =  -footing.depth; // depth under load bottom
	for (i=0; i<dividedLayers.length; i++) {
		midPoint = _top + 0.5*dividedLayers[i].thickness;
		if (midPoint > 0) {
		  compactionArray.push(cc*getIz(midPoint,footing, IzPeak, IzDepth, useModified)*dividedLayers[i].thickness/dividedLayers[i].Youngs);
		}
		else {
			compactionArray.push(0.0);
		}
		_top += dividedLayers[i].thickness;
	}
	
	// now compute settlements
	settlementArray = new Array(dividedLayers.length+1);
	_bottom = _top;
	_settlement = 0.0;
	// let's make settlement x and depth y for graph
	// make y negative
	// if point is above the footing
	settlementArray[dividedLayers.length] = [_settlement,-_bottom]; 
	for (i=dividedLayers.length-1; i>=0; i--) {
		  _settlement += compactionArray[i];
		  _bottom -= dividedLayers[i].thickness;
		  // hack to nullify floating point number slightly less than 0
		  _bottomTemp = -_bottom;
		  if (_bottomTemp > 0) _bottomTemp = 0.0;
		  settlementArray[i] = [_settlement,_bottomTemp];
	}
	
	return settlementArray;
	
	// integrate
//	_top =  -footing.depth; // depth under load bottom
//	totalSettlement = 0.0;
//	for (i=0; i<dividedLayers.length; i++) {
//		midPoint = _top + 0.5*dividedLayers[i].thickness;
//		if (midPoint > 0) {
//		  totalSettlement += getIz(midPoint,footing, IzPeak, IzDepth, useModified)*dividedLayers[i].thickness/dividedLayers[i].Youngs;
//		}
//		_top += dividedLayers[i].thickness;
//	}
//	totalSettlement *= C1*C2*delta_q;
	
//	return totalSettlement;
}