]> --Ndens=1e15 --Lcell=10.0e-2 --Temperature=1e-3 --Pwidth=0.4e-6 --delta1=0 --delta2=0 --delta3=0 --E1o=0 --E2o=1e2 --E3o=0 --E4o=0 realistic_Rb Eugeniy Mikhailov License GPL. Solving simplified Rb atom model with fields propagation along spatial axis Z with Doppler broadening. We assume four-wave mixing condition when w3-w4=w2-w1 i.e. fields E3 and E4 drive the same resonance as fields E2 and E1. * --------------- | F=1, 2P_3/2 > * \ \ * \ E3_r \ -------- | F=2, 2P_+1/2 > * \ E4_r \ / \ * \ \ / E2_l \ * \ / \ E1_l * | F=2, 2S_1/2 > -------------- \ * \ \ * \ \ * ------------- | F=1, 2S_1/2 > * We are solving dE/dz+(1/c)*dE/dt=i*eta*rho_ij, where j level is higher than i. Note that E is actually a Rabi frequency of electromagnetic field not the EM field in xmds terms it looks like dE_dz = i*eta*rhoij - 1/c*L[E], here we moved t dependence to Fourier space VERY IMPORTANT: all Rabi frequency should be given in [1/s], if you want to normalize it to something else look drho/dt equation. No need to renormalizes eta as long as its express through the upper level decay rate in the same units as Rabi frequency. &RbAtomConstansFile; section double v_thermal_averaged=0; // Maxwell distribution velocities range to take in account in [m/s] double V_maxwell_min = 0, V_maxwell_max = 0; // repopulation rate (atoms flying in/out the laser beam) in [1/s] const double gt=0.01 *(2*M_PI*1e6); // Larmor frequency double WL=0; // inner use variables double probability_v; // will be used as p(v) in Maxwell distribution ]]> 0 to provide range for Maxwell velocity distribution\n"); v_thermal_averaged=sqrt(k_boltzmann*Temperature/mass); // Maxwell distribution velocities range to take in account in [m/s] // there is almost zero probability for higher velocity p(4*v_av) = 3.3e-04 * p(0) V_maxwell_min = -4*v_thermal_averaged; V_maxwell_max = -V_maxwell_min; ]]> z E1 E2 E3 E4 probability_v probability_v_norm Maxwell_distribution_probabilities E1a E2a E3a E4a E_field Maxwell_distribution_probabilities Maxwell_distribution_probabilities_norm &RbChosenRhoFile; // this is because if r11=1, everything else is 0 and then every small increment // seems to be huge and adaptive solver makes smaller and smaller steps. // As quick and dirty fix I reshuffle initial population // so some of the population sits at the second ground level |2> // TODO: Fix above. Make the equation of motion for r11 // and express other level, let's say r44 // through population normalization ]]> &RbInitsFile; 100 density_matrix E_field_avgd &RbEquationsFile; Lt E_field density_matrix

E_field_avgd I1_out I2_out I3_out I4_out