Integración de Monte Carlo en MPI

Sea 4 dimensiones: x, y, z, t. Hay un toro de tamaño variable en el tiempo:

(R(t) - sqrt(x²+y²) )² +z² <= r(t)²

 
R(t): distancia del centro del tubo del toro al centro del toro (0,0,0)
r(t): radio del tubo

Dentro del toro se genera energia, con una potencia de P(x,y,z,t) Watts (=Joules/s)
Diga cuanta energia se genera desde el tiempo 0 al tiempo 10 para los siguientes R,r,P:

R(t)=100-t; r(t)=(t+1)²/2; P(x,y,z,t)=abs(z)  

#include <stdio.h >
#include <stdlib.h >
#include <time.h >
#include <math.h >
#include "mpi.h"

double R(double t){
  return 100.0 - 1.0*t;
}

double r(double t){
  return ((t+1)*(t+1))/2.0;
}

double P(double x, double y, double z, double t){
  return (double)abs(z);
}

int main(int argc, char *argv[])
{
  int nprocs, myproc;
  int i, n;
  double acum, sum_acum;
  double E, w, x, y, z, t, xx, yy, R_prima, r_prima, W;
  double t_i, t_f;
  double t0, t1;
  double x_max, x_min;
  double y_max, y_min;
  double z_max, z_min;
  double avg;

  /* Init MPI */
  MPI_Init( &argc, &argv );
  MPI_Comm_size( MPI_COMM_WORLD, &nprocs );
  MPI_Comm_rank( MPI_COMM_WORLD, &myproc );

  if( argc == 2 ) {      
    sscanf( argv[1], "%d", &n );
  }
  else {
    if (myproc == 0)
    {
      printf("Input the number of random points: (0 quit) ");fflush(stdout);
      scanf("%d",&n);     
    }
 }

  if( myproc == 0 ) printf("Calculating E using %d points\n", n);

  MPI_Bcast( &n, 2, MPI_INT, 0, MPI_COMM_WORLD);
  MPI_Barrier( MPI_COMM_WORLD );
  t0 = MPI_Wtime();

  E = 0.0;
  sum_acum = 0;
  acum = 0;
  t_i = 0.0;
  t_f = 10.0;
  x_min = y_min = -150.5; //-(R(t) + r(t))
  x_max = y_max = 150.5;  //R(t) + r(t)
  z_min = -60.5;  // -r(t)
  z_max = 60.5;   // r(t)
  srand( (unsigned)time( NULL )+myproc );
  
  for( i=myproc; i<n; i+=nprocs ) {
    t = ((double) rand()/RAND_MAX)*(t_f - t_i);    
    x = (((double) rand()/RAND_MAX)* (x_max - x_min)) - x_max;
    y = (((double) rand()/RAND_MAX)* (y_max - y_min)) - y_max;
    z = (((double) rand()/RAND_MAX)* (z_max - z_min)) - z_max;
    r_prima = r(t);
    R_prima = R(t);
    xx = ((double) R_prima - (double)sqrt(x*x + y*y));
    if(xx*xx + z*z <= r_prima*r_prima)
      acum += P(x, y, z, t); 
  }
  MPI_Allreduce( &acum, &sum_acum, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD );
  avg = (double)sum_acum / (double)n; 
  
  E = (t_f - t_i) * (x_max - x_min) * (y_max - y_min) * (z_max - z_min) * avg;

  t1 = MPI_Wtime();
  if( myproc == 0 ) {
    printf("E =  %f Joules\n", E);
    printf("Time = %f s.\n", t1-t0);
  }
  MPI_Finalize();
  return 0;
}

Message Passing Interface (MPI)

• M P I = Message Passing Interface
• MPI es una especificación para los desarrolladores y usuarios de las librerias de paso de mensajes . Por sí mismo, no es una libreria - sino más bien la especificación de lo que debe ser dicha libreria.
• MPI adopta  principalmente el modelo de programación paralela de paso de mensajes: los datos se mueven desde el espacio de direcciones de un proceso a la de otro proceso a través de las operaciones de cooperación en cada proceso.
• En pocas palabras, el objetivo de la MPI es proporcionar un estándar ampliamente utilizado para la escritura de programas de paso de mensajes. La interfaz intenta ser :
• práctico
• portátil
• eficiente
• flexible
• El estándar MPI ha pasado por una serie de revisiones, con la versión más reciente de las cuales MPI- 3. 

Modelo de programación:

• Originalmente, MPI ha sido diseñado para las arquitecturas distribuidas de memoria, que se estaban volviendo cada vez más popular en ese momento (1980 - principios de 1990).

 

• Como las tendencias de arquitectura cambian, las SMP de memoria compartida se combinaron a través de redes que crean memorias sistemas de memoria distribuida híbridos/compartidos.
• Ejecutores MPI adaptaron sus librerias para manejar ambos tipos de arquitecturas de memoria asociados a la perfección. También adaptado/desarrollaron maneras de manejar diferentes interconexiones y protocolos.

• Hoy, MPI funciona en prácticamente cualquier plataforma de hardware:

• memoria distribuida 
• memoria compartida 
• híbrido

Razones para el uso de MPI:

• Normalización - MPI es la única biblioteca de paso de mensajes que puede ser considerado un estándar. Cuenta con el apoyo de prácticamente todas las plataformas HPC. Prácticamente, se ha sustituido todas las bibliotecas de paso de mensajes anteriores.

• Portabilidad - Hay poca o ninguna necesidad de modificar el código fuente al puerto de la aplicación para una plataforma diferente que apoya (y es compatible con) el estándar MPI.

• Oportunidades de Rendimiento - implementaciones de proveedores deben ser capaces de explotar las características de hardware nativo para optimizar el rendimiento.

• Funcionalidad - Hay más de 440 rutinas definidas en MPI-3, que incluye la mayoría de los de MPI-2 y MPI-1.

• Disponibilidad - Una variedad de implementaciones están disponibles, tanto de proveedores y de dominio público.

Documentation:

• Documentación para todas las versiones del estándar MPI está disponible en:http://www.mpi-forum.org/docs/.

Integración de Monte Carlo

• El método de "Monte Carlo" es una técnica numérica que usa números pseudoaleatorios.
• La Integración de Monte Carlo estima el valor de una integral.
• Toma valores de la función en puntos aleatorios.
• El área (volumen) por el valor promedio estima la integral.

• ¿Cuándo se usa la Integración de Monte Carlo?
• Cuando la función a evaluar es complicada.
• Normalmente para integrales múltiples.
• La integración de Monte Carlo siempre funciona, pero toma mucho tiempo de CPU para obtener una respuesta más exacta.
• Pero, CPU son baratas.

Es decir, los métodos de integración de Montecarlo son algoritmos para encontrar una evaluación aproximada de una integral definida, normalmente de integrales múltiples. Los algoritmos deterministas de integración numérica, para aproximar la integral, evalúan la función en un conjunto de puntos correspondientes a una parrilla regular o en un conjunto de puntos predefinidos. 

En cambio, los métodos de Montecarlo eligen de forma aleatoria los puntos en los que se evaluará la función. La integración de Montecarlo forma parte de una familia de algoritmos llamados genéricamente métodos de Montecarlo.

Método:

• Recordemos el primer año de cálculo: El teorema del valor medio.

 

• La integral de la función f, sobre el volumen V, es igual a la media de f por V.
•  La integración de Monte Carlo esta dado por:
• Estimación del valor medio de la función.
• Multiplicado por el volumen. 

•  Estimación del valor medio de la función: 
• La media de una función es calculado por la media simple (también llamado promedio)

•  La varianza de la media simple es

• El estimado de la media de la función donde generamos un gran número de valores aleatorios de xi, y los usamos para calcular la media simple de la función.

 

•  Calculando la Integral de Monte Carlo

 

• EJEMPLO: Cálculo del número PI

/* Autor: Jason Martinez
 * Calculating number PI, with Monte Carlo Integation
 *            x^2 + y^2 <= 1
 * A = PI*r^2
*/

#include <stdio.h >
#include <stdlib.h >
#include <time.h >
#include <math.h >

const double pi = 3.14159265358979323846264338327950288419716939937510;

int main()
{
  double r, x, y;
  double avg;
  double sum, A;
  double x_max, x_min, y_max, y_min;
  int n;  //number of iterations
  int i;
  double pi_calc, err;

  scanf("%d", &n);

  sum = 0;
  A = 0;
  r = 1.0;
  x_max = y_max = 1.0;
  x_min = y_min = -1.0;
  
  for(i=0; i<n; ++i)
  {
    x = ((double)rand()/RAND_MAX)*(x_max - x_min) + x_min;
    y = ((double)rand()/RAND_MAX)*(y_max - y_min) + y_min;
    if(x*x + y*y <= r)
      sum += 1.0;
  }

  avg = sum / n; 

  A = (x_max - x_min)*(y_max - y_min) * avg;
  
  pi_calc = A/(r*r);
  printf("pi = %lf\n", pi_calc);
  err = fabs(pi - pi_calc);
  printf("error = %lf\n", err);
  return 1;
}

Introducción

Desde hace unos años las CPUs tienen ciertas limitaciones de hardware que hacen necesaria una refrigeración especial. En este punto la carrera para reducir el tamaño del silicio e incrementar la frecuencia del reloj ha terminado. En lugar de gastar grandes cantidades de silicio en el CPU para algoritmos avanzados que mejoren las instrucciones pre-fetch, CPUs más pequeños y más simples  son usados y hay espacio para más de una CPU en una misma oblea de silicio. Tenemos la CPU multi-núcleo que significa prácticamente varias CPUs en el mismo equipo.

Al principio, los núcleos de una CPU multi-núcleo eran más sencillos que el de un solo núcleo. Estos núcleos también operan en una frecuencia mucho más baja lo que significa que una aplicación diseñada para una operación única de la tarea tendría un enorme impacto en el rendimiento al pasar a un equipo nuevo, por primera vez.

La computación paralela se ha convertido en la corriente principal. Empezamos con una larga serie de conferencias sobre la computación paralela. Parecía que la gente quería saber sobre este tema, pero sobrevino que la computación paralela asusto a la gente. Hay una enorme brecha antes de que puedas llegar a ser un buen programador paralelo. Así como lo es para la programación orientada a objetos. Esto significa que los jefes de equipo y los arquitectos se encontraban en el mismo nivel que los programadores principiantes.

Añade a esto el hecho de que hay grandes cantidades de código ya escrito para una CPU de un solo núcleo y las ventajas que pueden obtenerse después de que al menos alguien lo re-escriban. Pero la más importante razón de peso para rechazar la computación paralela fue que es más fácil y más barato comprar otra máquina, que hacer el mejor de los núcleos del CPU. Esta fue una realidad que impulso la computación de la nube (Cloud Computing).

• Computador paralelo: Capaz de ejecutar varias instrucciones simultáneamente.

• Computación Paralela: Uso de varios procesadores trabajando juntos para resolver una tarea común:
• Cada procesador trabaja en una porción del problema.
• Los procesos pueden intercambiar datos, a través de la direcciones de memoria compartidas o mediante una red de interconexión

• Programación Paralela: Considera aspectos conceptuales y las particularidades físicas de la computación paralela.
• Objetivo: Mejorar las prestaciones mediante un buen aprovechamiento de la ejecución simultánea