Added the computation of Compton scattering using the Klein-Nishina cross-section (6b75e7e4) · Commits · polarization / slab

Makefile

+11 −2

Original line number	Diff line number	Diff line
		LDFLAGS = -L ${LIB_CFITSIO} -L ${LIB_GSL} -L ${LIB_MPI} -lcfitsio -lm -lgsl -lgslcblas -lmpi
		INCLUDE = -I ${PWD} -I ${HEADERS_CFITSIO} -I ${HEADERS_GSL} -I ${HEADERS_MPI} -I ${MC_INSTDIR}/../../include/
		LDFLAGS = -L ${LIB_GSL} -L ${LIB_CFITSIO} -L ${LIB_MPI} -lcfitsio -lm -lgsl -lgslcblas -lmpi
		INCLUDE = -I ${PWD} -I ${HEADERS_GSL} -I ${HEADERS_CFITSIO} -I ${HEADERS_MPI} -I ${MC_INSTDIR}/../../include/

		CFLAGS = -g -O3 -Wall ${INCLUDE}

		@@ -54,6 +54,15 @@ endif
		${CC} $^ -o $@ ${LDFLAGS}


		debug: CFLAGS += -pg
		debug: LDFLAGS += -pg
		debug: clean
		debug: start_iteration
		@echo ""
		@echo "✔ Compilazione completata in modalità DEBUG con gprof."
		@echo "Esegui ./start_iteration normalmente, poi usa:"
		@echo " gprof start_iteration gmon.out > profile.txt"
		@echo ""

		clean:
		rm -f main_program.o ${LOCAL_DIR_OBJ} ${MC_GENERAL_OBJ} ${MC_POLARIZATION_OBJ} ${MC_ELE_DISTRIBUTION_OBJ} ${MC_COORD_OBJ} ${MC_SEED_OBJ} ${MC_HYDRO_OBJ} ${MC_SLAB_OBJ} ${MC_SYNC_OBJ}

compute_results.c

+9 −10

Original line number	Diff line number	Diff line
		@@ -84,7 +84,6 @@ void compute_results(int method,
		fprintf(dat, "!Col2: Pdeg\n");
		fprintf(dat, "!Col3: I(u)/I(u=1)\n");


		if (disktau < 4)
		{
		for (jj = 0; jj < Nstep_mu; jj++)

compute_stokes.c

+6 −8

Original line number	Diff line number	Diff line
		@@ -36,13 +36,11 @@ void compute_stokes(double ene, double* k_lab, double* polvect_cart, stokes_para
		{
		if ((cos_theta <= array_ubounds[jj]) && (cos_theta >= array_ubounds[jj + 1]))
		{

		detector_axis(k_lab, Qp_cart, Qm_cart, Up_cart, Um_cart);

		Qs = pow(dot_prod(Qp_cart, polvect_cart), 2) - pow(dot_prod(Qm_cart, polvect_cart), 2);
		Us = pow(dot_prod(Up_cart, polvect_cart), 2) - pow(dot_prod(Um_cart, polvect_cart), 2);


		if (isnan(Qs))
		{
		printf("Error: Qs is Nan\n");

legendre_integration.c

+4 −6

Original line number	Diff line number	Diff line
		@@ -12,7 +12,6 @@ Gauss-Legendre quadrature
		/* [2(1-u^2)(1-u'^2) + u^2 u'^2] I_l(u) du' */
		/========================================================/


		double legendre_integration_A(gsl_spline2d* Spline_I2d, double tau, double u, double* u_prime, double* weights_u)
		{
		int status;
		@@ -38,7 +37,6 @@ double legendre_integration_A(gsl_spline2d* Spline_I2d, double tau, double u, do

		// printf("u %lf tau %lf Ik_l %5.4e\n", u_prime[ii], tau, Ik_l);


		if (Ik_l < -tolerance)
		{
		//("Warning: interpolated Ik_l < 0 for tau=%5.4f u_prime=%5.4f\n", tau, u_prime[ii]);

main_program.c

+55 −78

Original line number	Diff line number	Diff line
		#include "variables.h"
		#include "globals.h"
		#include "variables.h"
		#include <mpi.h>
		#include <stdio.h>


		char* set_filename(char* root, char* tau_char, char* seed_char, char* albedo, char* method);
		char* set_mc_filename(const char* root, double ktbb, double kte, double tau, int seed, double albedo,
		const char* method);
		char* set_mc_filename(const char* root, double ktbb, double kte, double tau, int seed, double albedo, const char* method);
		void smart_double_to_string(char* buffer, size_t size, double value);

		void int_to_string(char* buffer, size_t size, int value);



		int seed;
		double tau_c;
		double albedobase;
		@@ -20,8 +16,6 @@ double albedobase;
		char* polarization_file;
		char* diffusion_file;



		int kiter;
		int DiskPolarizationFlag;
		int Albedoflag;
		@@ -178,11 +172,8 @@ int main(int argc, char* argv[])
		{
		nstepangles = (int)numpar(argv[t]);
		}

		}



		/=========================================================/
		/Default values/
		/=========================================================/
		@@ -218,7 +209,6 @@ int main(int argc, char* argv[])
		nstepangles = 10;
		}


		/=========================================================/

		if (!disktau \|\| disktau == 0)
		@@ -289,7 +279,8 @@ int main(int argc, char* argv[])

		if (nphtot == 0)
		{
		printf("Please select the number of photons for MC simulation with parametere nphtot=nphtot\n");
		printf(
		"Please select the number of photons for MC simulation with parametere nphtot=nphtot\n");
		exit(1);
		}

		@@ -386,10 +377,6 @@ char* set_filename(char* root, char* tau_char, char* seed_char, char* albedo, ch
		return filename;
		}





		/======================================================================================/

		// Conversione intelligente dei double (senza decimali inutili)
		@@ -407,11 +394,7 @@ void int_to_string(char* buffer, size_t size, int value)
		snprintf(buffer, size, "%d", value);
		}




		char* set_mc_filename(const char* root, double ktbb, double kte, double tau, int seed, double albedo,
		const char* method)
		char* set_mc_filename(const char* root, double ktbb, double kte, double tau, int seed, double albedo, const char* method)
		{
		// Buffer per le conversioni
		char ktbb_str[20], kte_str[20], tau_str[20], albedo_str[20], seed_str[20];
		@@ -423,12 +406,9 @@ char* set_mc_filename(const char* root, double ktbb, double kte, double tau, int
		int_to_string(seed_str, sizeof(seed_str), seed);

		// Calcolo dimensione necessaria per il nome file
		size_t len = strlen(root)
		+ strlen("_ktbb") + strlen(ktbb_str)
		+ strlen("_kte") + strlen(kte_str)
		+ strlen("_tau")+ strlen(tau_str)
		+ strlen("_seed") + strlen(seed_str)
		+ strlen(".qdp") + 1;
		size_t len = strlen(root) + strlen("_ktbb") + strlen(ktbb_str) + strlen("_kte") +
		strlen(kte_str) + strlen("_tau") + strlen(tau_str) + strlen("_seed") +
		strlen(seed_str) + strlen(".qdp") + 1;

		if (strcmp(method, "mc") == 0)
		{
		@@ -464,6 +444,3 @@ char* set_mc_filename(const char* root, double ktbb, double kte, double tau, int

		return filename; // Ricorda di liberare con free()
		}