fibercryst example working : 1W -> 1.22W with 0.88W pumping

program2D/main.cc

@@ -174,11 +174,12 @@ void amplification(double dz,
     double c = 3e11;// mm / s
     //c = 3e8;// m / s
     double planck = 6.626e-34;
-    double emissionCrosssection = 7.7* 1e-18;//mm²
-//    emissionCrosssection = 7.7* 1e-20;//cm²
-//    emissionCrosssection = 7.7* 1e-24;//m²
+    double emissionCrosssection = 2.5* 1e-17;//mm²
     double upperLevelLifetime = 230e-6;
-    double Ntot = 1.3e+17; //( 1 / mm^3)
+    double Ntot = 1.39e+17; //( 1 / mm^3)
+    double lambdaPump = 808e-6;
+          // lambda   = 1064e-6;
+    //doping 0.1&
     //Ntot = 1.3e+26; //( 1 / m^3)
 
     //pumppower : W / mm³
@@ -187,36 +188,36 @@ void amplification(double dz,
     VariableFFT inversion(temp);
     VariableFFT pumpPhotons(temp);
     Function2d1<double, std::complex<double> > absolute(ABS);
-    // photonDensity : Intensity / energie_photon / dz , unit: 1 / (s mm^3), dz richtig hier?
-    photonDensity = absolute(varIn) * absolute(varIn)/ (planck * c / lambda) / dz;
+    // photonDensity : Intensity / energie_photon , unit: 1 / (s mm^2), dz richtig hier?
+    photonDensity = absolute(varIn) * absolute(varIn)/ (planck * c / lambda);
     // pumpPhotons : pumppower / energie_photon , unit: 1 / (s mm^3)
-    pumpPhotons = pumppower / (planck * c / lambda); // N / s mm³
+    pumpPhotons = pumppower / (planck * c / lambdaPump); // N / s mm³
 
-    std::cout << "pumppower "<< L_infty(pumppower) << std::endl;
-    std::cout << "inversion "<< L_infty(temp) << std::endl;
-    std::cout << "photonDensity "<< L_infty(photonDensity) << std::endl;
+//    std::cout << "pumppower "<< L_infty(pumppower) << std::endl;
+//    std::cout << "inversion "<< L_infty(temp) << std::endl;
+//    std::cout << "photonDensity "<< L_infty(photonDensity) << std::endl;
     // exp complex: normale exp funktion, die aber den realteil des arguments in den exponenten packe, also exp(real(arg)) und einen complexen wert (imag = 0) zurückgibt
     //noetig, weil multipliziert mit komplexer variable
     Function2d1<std::complex<double>,std::complex<double>> Exp(expComplex);
 
     //inversion : stimmt das so?
-    inversion = pumpPhotons / (emissionCrosssection * pumpPhotons * photonDensity * c + 1.0 / upperLevelLifetime + pumpPhotons / Ntot);
+    inversion = pumpPhotons / (emissionCrosssection * photonDensity  + 1.0 / upperLevelLifetime + pumpPhotons / Ntot);
 
-    std::cout << "pumppower "<< L_infty(pumppower) << std::endl;
-    std::cout << "photonDensity "<< L_infty(photonDensity) << std::endl;
-    std::cout << "inversion "<< L_infty(inversion) << std::endl;
-    std::cout << "pumpPhotons "<< L_infty(pumpPhotons) << std::endl;
+//    std::cout << "pumppower "<< L_infty(pumppower) << std::endl;
+//    std::cout << "photonDensity "<< L_infty(photonDensity) << std::endl;
+//    std::cout << "inversion "<< L_infty(inversion) << std::endl;
+//    std::cout << "pumpPhotons "<< L_infty(pumpPhotons) << std::endl;
     //temp : argument, welches in den exponenten zur verstärkung kommt
-    temp= emissionCrosssection * inversion * photonDensity * dz   ;
-    std::cout << "arg exp temp "<< L_infty(temp) << std::endl;
+    temp= emissionCrosssection * inversion * dz   ;
+    //std::cout << "arg exp temp "<< L_infty(temp) << std::endl;
     // e^(verstärkung)
     temp = Exp(temp);
-    std::cout << "exp temp "<< L_infty(temp) << std::endl;
+//    std::cout << "exp temp "<< L_infty(temp) << std::endl;
 
-        std::ofstream DATEIC;
-        DATEIC.open("varGAIN.vtk");
-        temp.Print_VTK(DATEIC);
-        DATEIC.close();
+//        std::ofstream DATEIC;
+//        DATEIC.open("varGAIN.vtk");
+//        temp.Print_VTK(DATEIC);
+//        DATEIC.close();
     //verstärkungsschritt. hier eventuell sqrt(temp), da es sich hier nicht um die intensität, sondern sqrt(intensität) handelt?
     varIn = varIn * temp;
 
@@ -396,11 +397,12 @@ void estimateBeamWaist(VariableFFT& varIn,
    Function2d1<double, std::complex<double> > absolute(ABS);
    Intensity = absolute(varIn) * absolute(varIn);
 
+   double powerHelper = product(absolute(Intensity),absolute(unity));
 
-    double power = product(Intensity,unity);
-    double medianX = product(Intensity * X * X,unity) / power;
-    double medianY = product(Intensity * Y * Y,unity) / power;
-    std::cout << "power = " << power <<std::endl;
+    double medianX = product(Intensity * X * X,unity) / powerHelper;
+    double medianY = product(Intensity * Y * Y,unity) / powerHelper;
+    powerHelper = powerHelper *varIn.getHx() * varIn.getHy();
+    std::cout << "power = " << powerHelper <<std::endl;
     std::cout << "beamwaistX = " << 2.0 * sqrt(medianX) <<std::endl;
     std::cout << "beamwaistY = " << 2.0 * sqrt(medianY) <<std::endl;
 }
@@ -635,7 +637,7 @@ void estimateBeamQuality(VariableFFT& varIn,
 int main(int argc, char** argv) {
     std::ofstream DATEI;
 
-    int n = 8;
+    int n = 6;
 
     double R1 = 100;
     double R2 = 100;
@@ -726,6 +728,8 @@ int main(int argc, char** argv) {
     deltaN.interpolate(vardDeltaNComplex,0.0);
 
     VariableFFT varIn(varE);
+    VariableFFT unit(varE);
+    unit = 1.0;
     VariableFFT varIntensity(varE);
     VariableFFT varPhase(varE);
     VariableFFT varTempX(varE);
@@ -865,6 +869,14 @@ int main(int argc, char** argv) {
 //                     varIn,
 //                     varE,Kx,Ky,X,Y,temp);
 
+
+    double totalPumpPower = 0;
+    Function2d1<double, std::complex<double> > absolute2(ABS);
+
+    temp = absolute2(varE)*absolute2(varE);
+    double powerHelper = product(absolute2(temp),absolute2(unit));
+    powerHelper = powerHelper *varE.getHx() * varE.getHy();
+    std::cout << "power before : "<<powerHelper << std::endl;
     for (int iter = 0 ; iter < zIterator.getSize();iter++)
     {
 
@@ -879,6 +891,7 @@ int main(int argc, char** argv) {
          pumppower.interpolate(varPumppowerComplex,0.0);
 
 
+         totalPumpPower += product(absolute2(pumppower),absolute(unit)) * varE.getHx() * varE.getHy() * zIterator.get_hz();
 
 
 //         deltaN = 0.0;
@@ -908,7 +921,7 @@ int main(int argc, char** argv) {
 
             std::ofstream DATEIA;
 
-            DATEIA.open("/local/er96apow/FFT_results/fibercryst_lowres_"+std::to_string(iter)+".vtk");
+            DATEIA.open("FFT_results/fibercryst_lowres_"+std::to_string(iter)+".vtk");
             varIntensity = absolute(varIn) * absolute(varIn);
             varIntensity.Print_VTK(DATEIA);
             DATEIA.close();
@@ -930,8 +943,12 @@ int main(int argc, char** argv) {
         zIterator.next();
     }
 
+    temp = absolute2(varIn)*absolute2(varIn);
+    powerHelper = product(absolute2(temp),absolute2(unit));
+        powerHelper = powerHelper *varE.getHx() * varE.getHy();
+        std::cout << "power after  : "<<powerHelper << std::endl;
 
-
+        std::cout << "total pump power : " << totalPumpPower << std::endl;
 //    VariableFFT varB(varA);
 /*
     Blockgrid2D* blockGridRec = varE.Give_blockgrid();