codice:
/* Introduciamo le strutture, ripassiamo l'I/O
/* e le funzioni, impariamo ad evitare gli effetti */
/* collaterali nell'uso di funzioni */
#include <stdio.h>
#include <math.h>
#define EXIT_SUCCESS 0
#define EXIT_FAILURE -9
#define EULER 0
#define EULER_CROMER 1
/* Notate che il punto e virgola alla fine */
/* delle prossime righe e' cruciale! */
struct phaseSpace {
double x;
double v;
};
/* I PROTOTIPI */
/* m_ nelle righe sotto ci ricorda che si */
/* tratta di variabili mute, li' solo per */
/* memoria (potremmo non metterle nel prototipo, */
/* ma indicare solo il tipo) */
/* in questo esempio usiamo p per i puntatori */
struct phaseSpace *initHarmonicOscillator(double m_x0,
double m_v0);
struct phaseSpace *euler(double m_dt, double m_omega2,
struct phaseSpace *m_pXAndV);
struct phaseSpace *eulerCromer(double m_dt, double m_omega2,
struct phaseSpace *m_pXAndV);
/* notate che energy() e' generica, */
/* mentre forceHarmonicOscillator() e' proprio per */
/* l'oscillatore armonico (da cui il nome) */
double energy(double m_omega2, struct phaseSpace *m_pXAndV);
double forceHarmonicOscillator(double m_omega2, double m_x);
double potentialEnergyHarmonicOscillator(double m_omega2,
double m_x);
double myReadDouble(char *m_printMessage);
long myReadLong(char *m_printMessage);
/*********************************/
int main(void) {
double x0, v0, dt, omega2, totalTime, energy0, energy1;
long numberOfSteps, algorithm, i;
struct phaseSpace xAndV;
/* I simboli # sono per far piacere a gnuplot, per */
/* esempio, o ad aiutarci ad usare un filtro per eliminare */
/* le righe che non ci serviranno nella grafica */
printf("# Integrazione Oscillatore Armonico\n");
printf("#v1.1.1, 30 settembre 2005\n");
/* queste myRead sono li' per evitarci un po' di codice... */
totalTime =
myReadDouble(
"Inserire il tempo totale di integrazione");
algorithm =
myReadLong("Inserire 0 per Eulero, 1 per Eulero-Cromer");
dt = myReadDouble("Inserire dt");
x0 = myReadDouble("Inserire x0");
v0 = myReadDouble("Inserire v0");
omega2 = myReadDouble("Inserire omega2");
/* Il cast lo facciamo noi esplicitamente, per tenere */
/* le cose sotto controllo!! */
numberOfSteps = (long)(totalTime/dt);
printf("# dt = %g tempo totale = %g\n",
dt, totalTime);
printf("# omega2 = %g numero passi = %d\n",
omega2, numberOfSteps);
/* una funzione senza effetti collaterali. quel */
/* che cambia lo cambia con = */
xAndV = *initHarmonicOscillator(x0, v0);
printf("# Condizioni iniziali x0 = %g v0 = %g\n",
xAndV.x, xAndV.v);
energy0 = energy(omega2, &xAndV);
printf("# L'energia al tempo t = 0 vale %g\n", energy0);
if (algorithm == EULER) {
printf("# Usiamo l'algoritmo di Eulero\n");
for (i = 0; i < numberOfSteps; i++) {
/* un'altra funzione senza effetti collaterali... */
xAndV = *euler(dt, omega2, &xAndV);
energy1 = energy(omega2, &xAndV);
printf("%g %g %g %g\n",
(double)i*dt, xAndV.x, xAndV.v, energy1 - energy0);
}
} else if (algorithm == EULER_CROMER) {
printf("# Usiamo l'algoritmo di Eulero-Cromer\n");
for (i = 0; i < numberOfSteps; i++) {
/* un'altra funzione senza effetti collaterali... */
xAndV = *eulerCromer(dt, omega2, &xAndV);
energy1 = energy(omega2, &xAndV);
printf("%g %g %g %g\n",
(double)i*dt, xAndV.x, xAndV.v, energy1 - energy0);
}
} else {
printf("# TUONI E FULMINI:");
printf("L'algoritmo n. %ld non e' implementato. Errore.\n",
algorithm);
exit(EXIT_FAILURE);
}
energy1 = energy(omega2, &xAndV);
printf("# L'energia al tempo t = %d vale %g\n",
numberOfSteps, energy1);
return EXIT_SUCCESS;
}
/*********************************/
struct phaseSpace *initHarmonicOscillator(double x0,
double v0) {
static struct phaseSpace xAndV;
xAndV.x = x0;
xAndV.v = v0;
return &xAndV;
}
/*********************************/
struct phaseSpace *euler(double dt, double omega2,
struct phaseSpace *pXAndVOld) {
static struct phaseSpace xAndVNew;
/* che succede se omettete le prime parentesi */
/* nella riga sotto ?*/
xAndVNew.v = (*pXAndVOld).v +
forceHarmonicOscillator(omega2, (*pXAndVOld).x) * dt;
/* perche' qui usiamo -> e non . ? */
xAndVNew.x = pXAndVOld->x + pXAndVOld->v * dt;
return &xAndVNew;
}
/*********************************/
struct phaseSpace *eulerCromer(double dt, double omega2,
struct phaseSpace *pXAndVOld) {
/* Qual e' la differenza fra Eulero ed Eulero-Cromer? */
static struct phaseSpace xAndVNew;
xAndVNew.v = (*pXAndVOld).v +
forceHarmonicOscillator(omega2, (*pXAndVOld).x) * dt;
xAndVNew.x = pXAndVOld->x + xAndVNew.v * dt;
return &xAndVNew;
}
/*********************************/
double energy(double omega2, struct phaseSpace *pXAndV) {
/* Calcoliamo in realta' 2 E / m */
double localEnergy;
localEnergy = (*pXAndV).v * (*pXAndV).v
+ potentialEnergyHarmonicOscillator(omega2, (*pXAndV).x);
return localEnergy;
}
/*********************************/
double forceHarmonicOscillator(double omega2, double x) {
return - omega2 * x;
}
/*********************************/
double potentialEnergyHarmonicOscillator(double omega2,
double x) {
return omega2 * x * x;
}
/*********************************/
long myReadLong(char *printMessage) {
long inputData;
printf("# %s\n",printMessage);
fflush(stdout);
scanf("%ld",&inputData);
return inputData;
}
/*********************************/
double myReadDouble(char *printMessage) {
double inputData;
printf("# %s\n",printMessage); fflush(stdout);
scanf("%lg", &inputData);
return inputData;
}
Rispondi quotando
Quindi sarebbe la stessa cosa se invece di utilizzare la funzione io inserissi direttamente nella main