## etendency1.m -- Compute the bulk energy tendency of the earth
## for constant albedo and greenhouse effect proportional to saturation
## water vapor column plus a constant (vapor plus co2).

function etendency1

  ## set constants
  fs = 1370.;
  albedo = 0.3;
  sigma = 5.67e-8;
  tcc = 5418.;     # tcc = L/R_v
  greenhouser = 33.;
  tr = 288.;

  ## plot prelims
  closeplot;
  hold off;
  grid;
  xlabel("ocean temperature (K)");
  ylabel("net flux imbalance (W/m^2)");
  title("Lines 1-6: 0% to 100% water vapor greenhouse at 288 K");
  plot([250,320],[0,0]);
  hold on;

  ## loop on greenhouse0
  for gpercent = [100:-20:0];

    ## compute tendency
    max = 69;
    oceantemp0 = 250;
    doceantemp = 1;
    oceantemp = oceantemp0 + doceantemp*[1:max];
    greenhouse0 = gpercent*greenhouser/100;
    greenhouse = greenhouse0 + (greenhouser - greenhouse0)*exp(tcc/tr - tcc ./ oceantemp);
    tendin = fs*(1 - albedo)/4;
    trad = oceantemp - greenhouse;
    for i = 1:max
      if (trad(i) < 0)
	trad(i) = 0;
      endif
    endfor
    tendout = sigma*trad.^4;
    tendency = tendin - tendout;

    ## plot results
    gset key left;
    plot(oceantemp,tendency);
    hold on;
  endfor

  ## do hard copy
  gset terminal postscript eps;
  gset output "etendency1.eps";
  replot;

endfunction
