diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 502af46..2ef4e44 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -17,6 +17,7 @@ jobs: matrix: compiler: ["gcc", "clang++"] generator: ["make", "ninja"] + fail-fast: false name: "Ubuntu compile.sh :: ${{ matrix.compiler }} :: ${{ matrix.generator }}" runs-on: ubuntu-latest steps: @@ -40,6 +41,7 @@ jobs: build_type: "Debug" exe_linker_flags: "-fsanitize=address" cxx_flags: "-fsanitize=address" + fail-fast: false name: "Ubuntu :: ${{ matrix.compiler }} :: ${{ matrix.build_type }} :: ${{ matrix.cxx_flags }} :: system ccfits cfitsio" runs-on: ubuntu-latest steps: diff --git a/manual.pdf b/manual.pdf index 3e9f4ea..d7b52aa 100644 Binary files a/manual.pdf and b/manual.pdf differ diff --git a/src/Dust.cpp b/src/Dust.cpp index 0b210b1..40e393e 100644 --- a/src/Dust.cpp +++ b/src/Dust.cpp @@ -2754,11 +2754,12 @@ void CDustComponent::preCalcMieScatteringProb() #pragma omp parallel for for(int a = 0; a < int(nr_of_dust_species); a++) { + // Init arrays of interp + avg_scattering_frac[a] = new interp[nr_of_wavelength]; + phase_pdf[a] = new interp[nr_of_wavelength]; + if(sizeIndexUsed(a)) { - // Init arrays of interp - avg_scattering_frac[a] = new interp[nr_of_wavelength]; - phase_pdf[a] = new interp[nr_of_wavelength]; for(uint w = 0; w < nr_of_wavelength; w++) { // Init pointer arrays @@ -4748,7 +4749,11 @@ void CDustComponent::getEscapePhotonMie(CGridBasic * grid, } else { - cout << "\nHINT: Photon package intensity or first scattering matrix element is zero!\n" << endl; + if(tmp_stokes.I() <= 0.0) + cout << "\nERROR: Photon package intensity is zero or negative!\n" << endl; + if(mat_sca(0, 0) <= 0.0) + cout << "\nERROR: First scattering matrix element is zero or negative!\n" << endl; + tmp_stokes.clear(); pp_escape->setStokesVector(tmp_stokes); return; @@ -4926,7 +4931,11 @@ void CDustComponent::miesca(photon_package * pp, uint a, CRandomGenerator * rand } else { - cout << "\nHINT: Photon package intensity or first scattering matrix element is zero!\n" << endl; + if(tmp_stokes.I() <= 0.0) + cout << "\nERROR: Photon package intensity is zero or negative!\n" << endl; + if(mat_sca(0, 0) <= 0.0) + cout << "\nERROR: First scattering matrix element is zero or negative!\n" << endl; + tmp_stokes.clear(); pp->setStokesVector(tmp_stokes); return; @@ -4958,7 +4967,7 @@ void CDustComponent::miesca(photon_package * pp, uint a, CRandomGenerator * rand run_counter++; } if(run_counter == 1000 || abs(root_phi) > 1e-10) - cout << "\nERROR: No phi found\n" << endl; + cout << "\nERROR: No scattering angle phi found!\n" << endl; phi = PI - gamma + phi; diff --git a/src/Dust.h b/src/Dust.h index cb8f214..d493e9d 100644 --- a/src/Dust.h +++ b/src/Dust.h @@ -124,7 +124,7 @@ class CDustComponent nr_of_mixtures = 0; calorimetry_type = 0; alignment = ALIG_PA; - phID = 1; + phID = 0; nr_of_dust_species = 0; nr_of_incident_angles = 0; nr_of_scat_theta = 0; diff --git a/src/Parameters.h b/src/Parameters.h index 1cb64ce..b1f465e 100644 --- a/src/Parameters.h +++ b/src/Parameters.h @@ -682,10 +682,11 @@ class parameters uint getPhaseFunctionID(uint i) const { + if(phIDs.size() == 1) + return phIDs[0]; if(i < phIDs.size()) return phIDs[i]; - else - return PH_ISO; + return PH_ISO; } double getFHighJ() const diff --git a/src/Synchrotron.cpp b/src/Synchrotron.cpp index f74f112..2185954 100644 --- a/src/Synchrotron.cpp +++ b/src/Synchrotron.cpp @@ -240,9 +240,15 @@ syn_param CSynchrotron::get_Power_Law_Parameter(double n_e, // additional correction for g_min>1 (see Reissl et al. 2018) double gmin_den = pow(g_min, 1 - p); +<<<<<<< HEAD res.j_I *= gmin_den; // res.j_I /= (g_min * g_min); // SMA: old was hardcoded for p = 3.0 res.j_Q *= gmin_den; // res.j_Q /= (g_min * g_min); // SMA: old was hardcoded for p = 3.0 res.j_V *= gmin_den; // res.j_V /= (g_min * g_min); // SMA: old was hardcoded for p = 3.0 +======= + res.j_I *= gmin_den; + res.j_Q *= gmin_den; + res.j_V *= gmin_den; +>>>>>>> tmp double du = sqrt(res.j_I * res.j_I - res.j_Q * res.j_Q); @@ -273,9 +279,9 @@ syn_param CSynchrotron::get_Power_Law_Parameter(double n_e, // * log(g_min)/tan_theta; // additional correction for g_min>1 (see Reissl et al. 2018) - res.alpha_I *= gmin_den; //res.alpha_I /= (g_min * g_min); // SMA: old was hardcoded for p = 3.0 - res.alpha_Q *= gmin_den; //res.alpha_Q /= (g_min * g_min); // SMA: old was hardcoded for p = 3.0 - res.alpha_V *= gmin_den; //res.alpha_V /= (g_min * g_min); // SMA: old was hardcoded for p = 3.0 + res.alpha_I *= gmin_den; + res.alpha_Q *= gmin_den; + res.alpha_V *= gmin_den; // converting back into SI; res.scale(); diff --git a/src/Synchrotron.h b/src/Synchrotron.h index 49d4b57..7f5a8f0 100644 --- a/src/Synchrotron.h +++ b/src/Synchrotron.h @@ -495,6 +495,7 @@ class CSynchrotron return Gamma((3. * p + 12.) / 12.) * Gamma((3. * p + 22.) / 12.) / (4. * (pow(g_min, 1. - p) - pow(g_max, 1. - p))); } + double getI_Q_p(double p) {