<?php   /* Code © Charles Chandler http://qdl.scs-inc.us/?top=18881 */   /* See http://qdl.scs-inc.us/?top=18943 for the description. */   function StefanBoltzmannLaw($emissivity, $meters2, $kelvins, $knownDepth = 1) { /* Stefan-Boltzmann Law -------------------- P = eoAT^4 where: P = total radiated power (watts) e = emissivity (=1 for BB radiators) o = Stefan's Constant A = radiating area (m^2) T = temperature (K) */ $stefansConstant = 5.670367 * pow(10, -8); /* Here's the correction for the effective depth of plasma BB radiators versus the classical BB surfaces ($knownDepth = 1). */ if (!$knownDepth) $knownDepth = pow($kelvins / 10000, 2.3219) * 4.905; return $emissivity * $stefansConstant * $meters2 * pow($kelvins, 4) * $knownDepth; }   function Acceleration($watts, $seconds, $kilograms) { /* Returns velocity of an object, in m/s, given the power used to accelerate it, the time, & the object's mass.   joules = .5 * kilograms * velocity^2 watts * seconds = .5 * kilograms * velocity^2 (watts * seconds) / (.5 * kilograms) = velocity^2 sqrt((watts * seconds) / (.5 * kilograms)) = velocity */ return sqrt(($watts * $seconds) / (.5 * $kilograms)); } // $mass_kg = $solSys->sun->mass_kg; $area_m2 = $solSys->sun->area_m2; $begSpeed = LU(11416, 'new_dp_velocity'); // imploding dusty plasma velocity: 258,961,135 m/s, or 86.38% c $nowSpeed = $begSpeed; $secsInMY = kSecsInYear * pow(10, 6); $imgFolder = '2ndParty/Images/Charles/LightCurve'; $steps = 10000; // number of iterations (max = 35000 before memory error) $sunAge = 380.6611474853; // from a previous run of the 'sun' $mode, and used by the 'earth' mode // // which numbers to run... // last one set wins $mode = '7500'; $mode = '22000'; $mode = '10000'; $mode = 'solarFlares'; $mode = 'sun'; $mode = 'earth'; $mode = 'earthFlares';   if ($mode == '22000') { $title = ''; $maxMegaYear = 5778; $begKelvins = 22000; $xPower = -3; $yPower = -3; } elseif ($mode == '10000') { $title = 'percent of total population'; $maxMegaYear = 5110; $begKelvins = 10000; $xPower = -3; $yPower = -3; } elseif ($mode == '7500') { $title = ''; $maxMegaYear = 9400; $begKelvins = 7500; $xPower = -3; $yPower = -3; } elseif ($mode == 'sun') { $title = ''; $maxMegaYear = 545; $begKelvins = 10000; $xPower = 0; $yPower = -3; $sunTemp = 5525; } elseif ($mode == 'solarFlares') { $title = 'percent of total population'; $maxMegaYear = 5110; $begKelvins = 10000; $xPower = -3; $yPower = -3; } elseif ($mode == 'earth') { $title = ''; $maxMegaYear = 600; $begKelvins = 10000; $mass_kg = $solSys->earth->mass_kg; $area_m2 = $solSys->earth->area_km2 * 1000 * 1000; $xPower = 0; $yPower = -3; } elseif ($mode == 'earthFlares') { $title = ''; $maxMegaYear = 600; $begKelvins = 10000; $mass_kg = $solSys->earth->mass_kg; $area_m2 = $solSys->earth->area_km2 * 1000 * 1000; $xPower = 0; $yPower = -3; } // $xLabel = ($xPower == -3) ? 'billions of years' : 'millions of years'; $yLabel = ($yPower == -3) ? 'thousands of kelvins' : 'kelvins'; // $relRootName = $imgFolder.ucfirst($mode); $relTextFile = $relRootName.'.txt'; if (file_exists($relTextFile)) { $pts = array(); foreach(explode(kLF, file_get_contents($relTextFile)) as $line) { if (trim($line)) { // skip the blank line at the end list($nowMY, $nowKelvins) = explode(kTAB, $line); $pts[] = Pt($nowMY, $nowKelvins); // The last-assigned value of $nowKelvins is used later. } } } else { $pts = array(Pt($xLabel, $yLabel)); // labels go in the first element $pts[] = Pt(0, $begKelvins); // add the initial temp at the 0 time step   // Use a logarithmic time step, with short steps at first, when the // temp is falling fast, and longer steps later, as it levels off. $logMaxSeconds = log($maxMegaYear * $secsInMY);   // If it's a star, 1/3 is nuclear power, // which doesn't diminish with power ouput. $nonNukePower = (BeginsWith($mode, 'earth')) ? 1 : (2 / 3);   // M A I N L O O P $nowKelvins = $begKelvins; $oldSeconds = 0; for ($i = 1; $i <= $steps; $i++) { $nowWatts = StefanBoltzmannLaw(1, $area_m2, $nowKelvins, false); $nowSeconds = exp(($i / $steps) * $logMaxSeconds); $timeSpan = $nowSeconds - $oldSeconds; $deceleration = Acceleration($nowWatts, $timeSpan, $mass_kg); $nowSpeed -= $deceleration * $nonNukePower; $nowKelvins = $begKelvins * ($nowSpeed / $begSpeed); $pts[] = Pt($nowSeconds / $secsInMY, $nowKelvins); $oldSeconds = $nowSeconds; // The last-assigned value of $nowKelvins is used later.   // Snag the age at which the Sun hits 5525 K. if (($mode == 'sun') and ($nowKelvins > $sunTemp)) $sunAge = $nowSeconds / $secsInMY; } // // Write the points to a text file. $str = ''; foreach($pts as $pt) $str .= $pt['x'].kTAB.$pt['y'].kLF; file_put_contents($relTextFile, $str); } // // Add the offset curves for the discrepancy between the Stefan-Boltzman // prediction and the actual temperatures in the K class due to flares. if (($mode == '10000') or ($mode == 'solarFlares')) { $pts1 = $pts2 = $pts3 = array(); $offY = 0; $begX = 650; $endX = 1200; $incX = ($endX - $begX) / 6; $incY = 850 / 6; $cnt = 0; foreach($pts as $pt) { $x = $pt['x']; $y = $pt['y']; if ($y > 5200) { $pts1[] = Pt($x, $y); } else { $pts4[] = Pt($x, $y); if (($x > $begX) and ($x < $endX)) { if ($x > $begX + ($incX * $cnt)) { $cnt++; $y += 4500; $offY += $incY; } } $pts3[] = Pt($x, $y - $offY); } if (($x > 1200) and ($y > $nowKelvins)) $pts2[] = Pt($x, $y - 850); } } elseif ($mode == 'earthFlares') { $pts1 = $pts2 = $pts3 = array(); $offY = 0; $begX = 8.6542077541337; $endX = 51.139340047215; $incX = ($endX - $begX) / 6; $incY = 192 / 6; $cnt = 0; foreach($pts as $pt) { $x = $pt['x']; $y = $pt['y']; if ($y > 5200) { $pts1[] = Pt($x, $y); } else { $pts4[] = Pt($x, $y); if (($x > $begX) and ($x < $endX)) { if ($x > $begX + ($incX * $cnt)) { $cnt++; $y += 2700; $offY += $incY; } } $pts3[] = Pt($x, $y - $offY); } if (($x > $endX) and ($y > $nowKelvins)) $pts2[] = Pt($x, $y - 192); } } // // GraphPoints is just a QDL class for generating a PNG file, plotting // the points, labeling the axes, and adding annotations as specified. // Without it, you can still plot the points using the text file (above). $plotWide = 640 * 2; $plotHigh = 480 * 2; $marginLeft = 52 * 2; $marginRight = 42 * 2; $marginUpper = 20 * 2; if ($title) $marginUpper = 90; $marginLower = 52 * 2; $gr = 180; // gray // $img = new GraphPoints($title, $plotWide, $plotHigh, 255, 255, 255, true, 24); $img->SetTickLengths(8, 16, 20); $img->SetLineThickness(2); $img->Margins($marginLeft, $marginRight, $marginUpper, $marginLower); $img->SetPowers($xPower, $yPower); if ($mode == '10000') { $img->xMin = 0; $img->xMax = $maxMegaYear * pow(10, $xPower); $img->yMin = 2400 * pow(10, $yPower); $img->yMax = $begKelvins * pow(10, $yPower); $img->AddPoints($pts4, true, 0, 0, 0); $img->AddPoints($pts1, true, 0, 0, 0); $img->AddPoints($pts2, true, 255, 0, 0); $img->AddText(2000, 4200, 'stefan-boltzmann prediction', 'left'); $img->AddText(1300, 2550, 'required for ~76.45% in M class', 'left', 255, 0, 0); } elseif ($mode == 'solarFlares') { $img->xMin = 0; $img->xMax = $maxMegaYear * pow(10, $xPower); $img->yMin = 2400 * pow(10, $yPower); $img->yMax = $begKelvins * pow(10, $yPower); $img->AddPoints($pts4, true, $gr, $gr, $gr); $img->AddPoints($pts1, true, 0, 0, 0); $img->AddPoints($pts3, true, 0, 0, 255); $img->AddPoints($pts2, true, 255, 0, 0); $img->AddText(2000, 4200, 'stefan-boltzmann prediction', 'left', $gr, $gr, $gr); $img->AddText(1300, 2550, 'after energy loss to flares', 'left', 255, 0, 0); $img->AddText( 710, 9600, 'flares', 'left', 0, 0, 255); } elseif ($mode == 'earthFlares') { $img->xMin = 0; $img->xMax = $maxMegaYear * pow(10, $xPower); $img->yMin = 1500 * pow(10, $yPower); $img->yMax = $begKelvins * pow(10, $yPower); $img->AddPoints($pts4, true, $gr, $gr, $gr); $img->AddPoints($pts1, true, 0, 0, 0); $img->AddPoints($pts3, true, 0, 0, 255); $img->AddPoints($pts2, true, 255, 0, 0); $img->AddText(150, 3100, 'stefan-boltzmann prediction', 'left', $gr, $gr, $gr); $img->AddText(150, 1900, 'after energy loss to flares', 'left', 255, 0, 0); $img->AddText( 10, 8200, 'flares', 'left', 0, 0, 255); } elseif ($mode == 'earth') { $img->xMin = 0; $img->xMax = $maxMegaYear * pow(10, $xPower); $img->yMin = 1500 * pow(10, $yPower); $img->yMax = $begKelvins * pow(10, $yPower); $img->AddPoints($pts, true, 0, 0, 0); $img->AddText(150, 3100, 'stefan-boltzmann prediction', 'left'); $img->AddLine($sunAge, 10000, $sunAge, $img->yMin / pow(10, $yPower)); } else { $img->AddPoints($pts, true, 0, 0, 0); } // // horizontal lines for stellar classes if ($mode == 'sun') { $img->AddLine(0, 7500, $maxMegaYear, 7500); $img->AddText($maxMegaYear * .97, ( 7500 + 10000) / 2, 'A'); $img->AddLine(0, 6000, $maxMegaYear, 6000); $img->AddText($maxMegaYear * .97, ( 6000 + 7500) / 2, 'F'); $img->AddLine(0, 5200, $maxMegaYear, 5200); $img->AddText($maxMegaYear * .97, ( 5200 + 6000) / 2, 'G'); $img->AddLine($sunAge - 11, 5525, $sunAge + 11, 5525, 1, 2); $img->AddLine($sunAge, $sunTemp - 120, $sunAge, $sunTemp + 120, 1, 2); echo 'sun age: '.$sunAge.'<br />'; } else { if (($mode != '10000') and ($mode != 'solarFlares') and (!BeginsWith($mode, 'earth'))) $img->AddLine(0, 10000, $maxMegaYear, 10000); $img->AddText($maxMegaYear * .97, (10000 + 30000) / 2, 'B'); $img->AddLine(0, 7500, $maxMegaYear, 7500); $img->AddText($maxMegaYear * .97, ( 7500 + 10000) / 2, 'A'); $img->AddLine(0, 6000, $maxMegaYear, 6000); $img->AddText($maxMegaYear * .97, ( 6000 + 7500) / 2, 'F'); $img->AddLine(0, 5200, $maxMegaYear, 5200); $img->AddText($maxMegaYear * .97, ( 5200 + 6000) / 2, 'G'); $img->AddLine(0, 3700, $maxMegaYear, 3700); $img->AddText($maxMegaYear * .97, ( 3700 + 5200) / 2, 'K'); $img->AddLine(0, 2400, $maxMegaYear, 2400); $img->AddText($maxMegaYear * .97, ( 2400 + 3700) / 2, 'M'); if ($mode == 'earthFlares') $img->AddLine($sunAge, 10000, $sunAge, $img->yMin / pow(10, $yPower)); } // // Hash marks for known star population percentages. if (($mode == '10000') or ($mode == 'solarFlares')) { $h = (($img->yMax - $img->yMin) / pow(10, $yPower)) * .015; $x1 = $maxMegaYear; $x2 = $x1 - ((76.45 / 100) * $maxMegaYear); $x3 = $x2 - ((12.10 / 100) * $maxMegaYear); $x4 = $x3 - (( 7.60 / 100) * $maxMegaYear); $x5 = $x4 - (( 3.00 / 100) * $maxMegaYear);   $above = 450; $img->AddLine( 0, 2400, 0, $begKelvins - $h + $above); $img->AddLine($x1, 2400, $x1, $begKelvins - $h + $above); $img->AddLine($x2, 2400, $x2, $begKelvins - $h + $above); $img->AddLine($x3, 2400, $x3, $begKelvins - $h + $above); $img->AddLine($x4, 2400, $x4, $begKelvins - $h + $above); $img->AddLine($x5, 2400, $x5, $begKelvins - $h + $above);   $img->AddText( .5, $begKelvins + $h + 400, '.6' ); // A $img->AddText(($x4 + $x5) / 2, $begKelvins + $h + 80, '3' ); // F $img->AddText(($x3 + $x4) / 2, $begKelvins + $h + 80, '7.6' ); // G $img->AddText(($x2 + $x3) / 2, $begKelvins + $h + 80, '12.1' ); // K $img->AddText(($x1 + $x2) / 2, $begKelvins + $h + 80, '76.45'); // M } // // white background $img->Indexed(true, 255); $img->Create($relRootName.'_wbg.png'); $img->DeIndex(); // // black background $img->Invert(); $img->Brighten(50); $img->Indexed(true, 255); $img->Create($relRootName.'_bbg.png'); $img->Destroy(); // // preview echo '<img src="'.$relRootName.'_bbg.png" width="640" height="480" />';   ?>