Ende

Mit diesem Post ist dieser Blog offiziell beendet. Die Swiss Ephemeris unterstützt seit der Version 2.03 das Multithreading.

Die Implementierungsidee war sehr einfach und kommt ganz ohne das in diesem Blog skizzierte grundsätzliche Refactoring aus: Globale Variablen, die pro Thread separat behandelt werden müssen, bekommen einfach bei der Deklaration den Storage class specifier _thread und gehören somit zur Speicherklasse "TLS" (für thread-local storage), der von C-Compilern schon seit gut anderthalb Jahrzehnten für die gängigen Zielarchitekturen unterstützt wird. Daten von dieser Speicherklasse werden bei der Erzeugung eines neuen Threads dupliziert, so dass der Thread über eine eigene Version von ihnen verfügt und sich beim Ändern der Inhalte nicht mit anderen Threads in die Quere kommt. In der Swiss Ephemeris war das im wesentlichen die globale Variable swed - eine Sammelstruktur, die als Zwischenspeicher für alle Art von Daten verwendet wird, die den einzelnen Aufruf einer Funktion der Swiss Ephemeris überleben sollen.

Wäre ich selbst mehr als nur freizeitmässig dem Thema zugeneigt gewesen, hätte mich schon 2013 die folgende Antwort des Users skylendar in der Swiss Ephemeris-Mailingliste nachdenklich stimmen müssen. Ich hatte gefragt, ob sich schon einmal andere Entwickler dieses Themas angenommen hätten und bekam von jenem skylendar die lapidare Antwort, er hätte schon vor langer Zeit eine mehrfädige Version der Swiss Ephemeris geschrieben, diese jedoch mangels Interesse in der Mailingliste nicht weiter verbreitet. Da ich damals nur den Lösungsweg gesehen hatte, das Programm total umzubauen, um die globalen Variablen aus dem Code zu eliminieren, fragte ich ihn, wie er dieses Kunststück vollbracht habe. Seine - wiederum lapidare - Antwort hätte mich damals nachdenklich stimmen müssen. Stattdessen überlas ich sie:

Delivered-To: ruediger.plantiko@gmail.com
Received: by 10.194.135.203 with SMTP id pu11csp79435wjb;
        Fri, 1 Mar 2013 08:04:14 -0800 (PST)
Message-ID: 
In-Reply-To: 
User-Agent: eGroups-EW/0.82
From: "skylendar" 
Sender: swisseph@yahoogroups.com
Mailing-List: list swisseph@yahoogroups.com; contact swisseph-owner@yahoogroups.com
Delivered-To: mailing list swisseph@yahoogroups.com
Date: Fri, 01 Mar 2013 16:04:07 -0000
Subject: [SWISSEPH] Re: Is libswe thread-safe?

--- In swisseph@yahoogroups.com, Rüdiger Plantiko  wrote:
>
> Skylendar,
>
> I missed this. Did you rewrite
> sweph.c?
> If not, how did you avoid the usage of global and static variables?
>
> Regards,
> Rüdiger

I used the thread keyword wherever it is necessary.

Da das Multithreading in der Mailingliste aber auch weiterhin immer wieder zum Thema wurde, postete skylendar am 4. September 2015 schliesslich doch seinen Patch, den er damals noch auf der Version 1.71 aufgesetzt hatte.

Es zeigte sich: wenn man diesen Weg verfolgt, bringt man sich zwar um eine an sich gute Gelegenheit, seinen Code zu refaktorisieren, aber man kann mit vergleichsweise wenigen Änderungen die Mehrfädigkeit erreichen.

Ich folgte seinem Weg für das damals aktuelle Release 2.02.01:

• Ich änderte alle globalen Variablen in const, auf die im Programm nur lesend zugegriffen wurde. Davon gab es eine ganze Reihe (etwa Koeffizienten oder sonstige Konstanten aus Reihenentwicklungen).
• Alle globalen Variablen, die nicht nur gelesen, sondern auch verändert wurden, bekamen das neue Attribut TLS (ein Makro, hinter dem sich das Schlüsselwort __thread verbarg.
• Dasselbe für die zahlreichen statischen Variablen, die sich ja nur durch die eingeschränkte Sichtbarkeit (File oder Function) von den globalen unterscheiden, aber derselben Verletzlichkeit bei mehrfädiger Programmausführung unterliegen wie globale Variablen.

Mit folgendem kleinen C99-Testprogramm, das die POSIX Schnittstelle <pthread.h> verwendet, konnte ich auf meinem Rechner mit 2*4 Prozessorkernen eine Laufzeitverbesserung um den Faktor 5 für die zufällige Berechnung von Planetenpositionen erreichen. Wie ich aber schon früher hier begründete, geht es nicht nur um die Laufzeitverbesserung durch effizientere Ressourcennutzung. Ein anderer Gewinn ist, dass mehrfädige Clientprogramme (und Plattformen wie z.B. mod_perl oder Java Server Pages) die Swiss Ephemeris nur dann sicher und ohne Extraprogrammierung zur Synchronisierung (etwa mit Mutexen) benutzen können, wenn diese selbst mehrfägig ist. Daher erweitert die Mehrfädigkeit auch die Einsatzmöglichkeiten der Swiss Ephemeris.

/*
 * Test program for multithreading
 * Performs many swe_calc() executions in several threads simultaneously
 *
 * Syntax: 'test <number of threads> [ <calcs per thread> [ <repetitions> [ <trace mode> ] ] ]'
 *
 */

#include <stdio.h>
#include <stdlib.h>
#include <pthread.h>
#include <sys/time.h>
#include <stdbool.h>
#include "swephexp.h"
#include "sweph.h"

typedef struct {
  double jd,l;
  int ipl;
  bool ok;
  } calc_data;

typedef struct {
  calc_data *m;
  int length;
} calc_data_array;

double get_time();
void *get_longitude( void* arg );
void *get_longitudes( void* arg );
void prepare_random_arguments( int nthreads, int calcs_per_thread, calc_data_array* data );
void parse_cmdline_args( int argc, char* argv[], int* nthreads, int *calcs_per_thread, int* repetitions, bool* trc);
void print_trace( int index,int nthreads, int calcs_per_thread, calc_data_array *data);


int main( int argc, char* argv[] ) {
  
  int nthreads         = 1;
  int calcs_per_thread = 1;
  int repetitions      = 100;
  bool trc             = false;
  double total_time    = 0;

// Inititialize pseudo random generator with system time as seed
  srand( time(NULL) );

  parse_cmdline_args( argc, argv, &nthreads, &calcs_per_thread, &repetitions, &trc);
 
// Allocate calc data array on stack
  calc_data _data[nthreads][calcs_per_thread];
  calc_data_array data[nthreads];
  for (int i=0;i<nthreads;i++) {
    data[i].m      = _data[i];
    data[i].length = calcs_per_thread;
  }

  pthread_t tid[nthreads];
 
  for (int i=0; i<repetitions; i++ ) {

// Prepare random argument chunks for swe_calc
    prepare_random_arguments( nthreads, calcs_per_thread, data );

// Compute a chunk of nthread*calcs_per_thread data
    double begin = get_time( );
    if (nthreads == 1) {
      // Direct call for comparison
      get_longitudes( &data[0] );
      }
    else {
      // ... or multithreaded
      for (int j=0;j<nthreads;j++) {
        pthread_create( &tid[j], NULL, get_longitudes, (void *) &data[j] );
      }
      for (int j=0;j<nthreads;j++) {
        pthread_join( tid[j], NULL);
      }
    }
    double end = get_time();
    total_time += (end - begin);

// Details, if wanted
    if (trc) print_trace(i,nthreads,calcs_per_thread,data);

  }  

// Execution times summary
  printf( 
      "%d times (%d threads, in chunks of %d), "
      "done in %lf msec total (=%lf /each) \n",
      repetitions,
      nthreads,
      calcs_per_thread,
      total_time*1000,
      total_time*1000/(repetitions*nthreads*calcs_per_thread)
      );

} 

// The "worker" function for a chunk of calculations, executed simultaneously in several threads
void *get_longitudes( void* arg) {
  calc_data_array *data = (calc_data_array*) arg;
  for (int i=0;i<data->length;i++) {
    get_longitude(&data->m[i]);
  }
  return NULL;
}

// The "worker" function for a single calculation
void *get_longitude( void* arg ) {
  calc_data *data = (calc_data*) arg;
  char serr[255];
  int  iflag = 0, iret = 0;
  double xx[6] = {};
  iret = swe_calc( data->jd, data->ipl, iflag, xx, serr );
  data->ok = (iret == iflag);
  data->l = xx[0];
  return NULL;
  }


// Use pseudo random generator to generate some Julian dates and planet numbers
void prepare_random_arguments( int nthreads, int calcs_per_thread, calc_data_array* data ) {
  for (int i=0; i<nthreads;i++)
    for (int j=0; j<calcs_per_thread;j++){
      data[i].m[j].jd  = J1900 + (rand() % 10000); // range: 10000 days from 1.1.1900 onwards
      data[i].m[j].ipl = rand( ) % 10;
    }
}


void parse_cmdline_args( int argc, char* argv[], int* nthreads, int *calcs_per_thread, int* repetitions, bool* trc) {

  if (argc == 1) {
    printf( "Syntax: 'test <number of threads> [ <calcs per thread> [ <repetitions> [ <trace mode> ] ] ]'\n");
    exit(4);
  }
  if (argc > 1) sscanf( argv[1], "%d", nthreads    );
  if (argc > 2) sscanf( argv[2], "%d", calcs_per_thread );
  if (argc > 3) sscanf( argv[3], "%d", repetitions );
  if (argc > 4) *trc = true; // *Any* fourth argument activates trace mode

}

double get_time() {
  struct timeval tv;
  gettimeofday(&tv,NULL);
  return tv.tv_sec + tv.tv_usec / 1000000.0 ;
  }

void print_trace( int index, int nthreads, int calcs_per_thread, calc_data_array *data) {
for (int k=0,kmax = nthreads*calcs_per_thread;k<kmax;k++) {
  int i = k / calcs_per_thread,
      j = k % calcs_per_thread;
  printf( "%d / %d / %d : %lf %d -> %lf %s\n",
          index,i,j,
          data[i].m[j].jd,
          data[i].m[j].ipl,
          data[i].m[j].l,
          data[i].m[j].ok ? "" : "(problem)");
  }
}