Pull to refresh

Comments 18

Что-то у меня под Windows собрать до работоспособного состояния пример extgcd.c, что там идёт в комплекте не получилось. Студии отчаянно не нравится заголовочник divll_pre.h, а gcc (minGW64) в общем собирает, и оно даже запускается, но крешится где-то в районе gp_read_stream(stdin). Я верю, что под Линуксом оно всё ОК, но вот Windows что-то не срастается.

Я пробовал и по-простому "gcc extgcd.c -llibpari -o extgcd" и так, как интернет рекомендует "gcc extgcd.c -Wall -Wextra -pedantic -llibpari -std=c11 -g -o extgcd", но нет.

Я так понял вы решили начать по хардкору, написав на Си с использованием libpari. Тут, действительно, потребуется больше времени.

Проще всего скомпилировать в библиотеку самую тяжёлую часть и загружать из pari. И даже это нужно делать только тогда, когда вся логика уже безукоризненно отлажена, а вычисления реально занимают часы или дни.

Я так понял вы решили начать по хардкору, написав на Си с использованием libpari. 

Это я просто "с конца" начал. Так то я в LabVIEW программирую, и хотелось бы (допустим) внешнюю математику подключить библиотекой (суть DLL), иначе pari останется "вещью в себе". Я, конечно, могу (вероятно) через командную строку с ней общаться, но хотелось бы именно библиотекой. Там ещё надо будет посмотреть как данные пробрасывать, но это дело техники - если я успешно соберу консольное приложение с libpari, то соберу и библиотеку. Там ещё могут настичь подводные камушки с параллельностью, так как LabVIEW "параллельна" сама по себе, и в некоторых редких случаях внешняя многопоточная библиотека может доставить головную боль, но и это решаемо. А "пристреливаться" и отлаживаться, я конечно не в библиотеке буду. У меня в общем практической задачки в руках нет, да и дополнительная зависимость - не всегда хорошо, тут просто чисто любознательный интерес для расширения кругозора.

CC = cc
INCDIR = /usr/include/pari
LIBDIR = /usr/lib64
CFLAGS = -O -I$(INCDIR) -L$(LIBDIR)
$(CC) $(CFLAGS) -o trans2c trans2c.c -lpari -lm

Я нашёл такую рекомендацию где-то в сети. Правда, для Linux.

А, не, там что-то с исходником не так. Я дома ещё раз проверил на совсем простом примере, всё ОК.

Вот код в test.c, я хз сколько там надо в pari_init() заказывать, надо доки почитать, пусть будет с запасом:

#include "pari/pari.h"

int main()
{
	pari_init(8000000,500000);
	pari_printf("Pari test - begin\n");

	GEN ver = pari_version(), Major = gel(ver,1), minor = gel(ver,2), p = gel(ver,3);
	pari_printf("Pari Version M = %Ps m = %Ps p = %Ps\n", Major, minor, p);
    printf("pari-%ld.%ld.%ld\n", itos(Major), itos(minor), itos(p));

	GEN M, D;
	M = mathilbert(5);
	D = det(M);
	pari_printf("M = %Ps\nD = %Ps\n", M, D);
	
	pari_close();
	printf("Pari test - end\n");	
	return 0;
}

Компилируем просто:

gcc test.c -l pari -o test

Вывод из test.exe:

Pari test - begin
Pari Version M = 2 m = 17 p = 0
M = [1, 1/2, 1/3, 1/4, 1/5; 
	1/2, 1/3, 1/4, 1/5, 1/6; 
	1/3, 1/4, 1/5, 1/6, 1/7; 
	1/4, 1/5, 1/6, 1/7, 1/8; 
	1/5, 1/6, 1/7, 1/8, 1/9]
D = 1/266716800000
pari-2.17.0
Pari test - end

И да, детерминант матрицы Гильберта 5x5 считает верно.

Параллельные вычисления теперь доступны каждому. С языком GP

Прочитал, но не понял это утверждение, оно написано так, словно параллельные вычисления были недоступны. В статье вы показали метод parfor который полностью эквивалентен Parallel.For/Parallel.Foreach в C# которые доступны всем давно. И так получается что писать ту самую параллельность вам придётся самому - а это очень специфичная тема и "познаниями бейсика" вы тут не обойдётесь.

Пока мне не удалось увидеть в PARI ничего интересного по сравнению с C#, и если последний будет легко изучать тем кто видел C/C++/Java, то PARI - это свой специфичный синтаксис.

В общем как-то статья не убедила, я, если честно, думал вы расскажите о чем-то, действительно уровня бейсика, где мне не придётся думать о памяти, общих переменных и т.п. и я просто отдам код, а оно само выполнится параллельно. Но это же не так.

Параллельные вычисления в языках общего назначения - это сложно и если нужно за полчаса написать программу, то вы не сделаете этого никогда не запускав C# до этого. В PARI достаточно, чтобы каждая процедура была автономной для параллельности.

А c PARI сделаете. А учитывая, что PARI поддерживает неограниченную точность и огромное количество математических функций, то задача ещё упрощается.

Параллельные вычисления в языках общего назначения - это сложно и если нужно за полчаса написать программу, то вы не сделаете этого никогда не запускав C# до этого. В PARI достаточно, чтобы каждая процедура была автономной для параллельности.

и в чем разница между parfor и Parallel.For/Parallel.Foreach в C#? Где конкретно в статье пишется что в PARI всё просто?

А c PARI сделаете

ну пока я вижу чтобы начать пользоваться PARI его нужно точно так же изучать, в тогда разница?

А учитывая, что PARI поддерживает неограниченную точность и огромное количество математических функций, то задача ещё упрощается.

Не слишком понял как это следует из статьи и как это обстоятельство, например неограниченная точность, влияет на простоту параллельного программирования.

Кстати, о самой фразе, которую вы процитировали. Это разные предложения.

Язык GP просто облегчает задачу, когда нужно что-то по-быстрому посчитать.

Я, например, могу написать параллельный код на Go или Python, но на GP для меня будет проще, если мне не нужны какие-то тонкости.

Кстати, о самой фразе, которую вы процитировали. Это разные предложения.

И?

Язык GP просто облегчает задачу, когда нужно что-то по-быстрому посчитать.

Каким образом?

Я, например, могу написать параллельный код на Go или Python, но на GP для меня будет проще, если мне не нужны какие-то тонкости.

Это может быть лично вашей спецификой, вопрос о конкретной статье и конкретных заявлениях автора - вас. Из статьи ничего подобного не следует.

Фактически, я часто просто вывожу в критической секции результаты параллельных функций в поток вывода или файл, а потом анализирую. Очень просто. Никаких заморочек с общими переменными. Вообще общие переменные использую только для чтения, как константы.

Но в критической секции можно результаты функций в какой-то глобальный массив записывать.

Фактически, я часто просто вывожу в критической секции результаты параллельных функций в поток вывода или файл, а потом анализирую. Очень просто

бла-бла-бла, это просто слова, вы статью же для чего-то писали, даже в конце вроде бы сослались на то - ради чего "Пробуйте и, может быть, вы разделите мой энтузиазм!".

У вас это энтузиазм, но статья не передаёт должной информации, чтобы этот энтузиазм появился у других. Вы точно так же как и для C# озвучили проблемы параллельных вычислений, точно такие же как и в C#, предложили (свели) всю параллельность к parfor - так что в PARI есть такое особенное? Я вам озвучил с какой надеждой я зашел читать статью, но я не увидел этого. Вы мне предложили запихать цикл в параллельный код, что я и сейчас могу делать в C#, причем очень много задач невозможно решить параллельным циклом не передавая информацию между итерациями. Это огромная проблема параллельных вычислений в целом и PARI ее не решает, вы же сами об этом пишете.

Никаких заморочек с общими переменными

И где это в статье? У вас вся параллельность, которая указана в заголовке, уместилась в один маленький раздел и в один пример.

Вообще общие переменные использую только для чтения, как константы.

Я точно так же это делаю в C#, но это сильное ограничение, например при обработке графов, когда у вас меняется вес в одной итерации и его нужно учитывать в другой - это огромная проблема параллельности вообще. Как её решает PARI?

Но в критической секции можно результаты функций в какой-то глобальный массив записывать.

Я могу это же делать в C#, собираю параллельную обработку в глобальный массив, а потом опять считаю, так как я не могу сразу работать с параллельными данными во время исполнения. Так в чем тут выигрыш PARI? Вы это как-то объяснили?

Повторю еще раз - parfor это точно такая же конструкция как и Parallel.For/Parallel.Foreach в C#, почему я должен бросить знакомый C#, изучать PARI? Вы же наверное собирались побудить людей к этому своей статьёй?

Я добавил раздел 0, где написано кому может понадобиться PARI. А также добавил ещё раздел 4.5, где обзорно пробежался по другим функциям параллельности в PARI.

Если сравнивать с C#, то порог вхождения несопоставим с PARI и значительно выше у C#.

C# - это язык общего назначения, а PARI - вещь заточенная для быстрого решения математических/вычислительных проблем.

спасибо, то что нужно.

прочитав раздел 0 в таком виде я просто закрыл бы статью и даже не писал бы комментарий, так как это не для меня.

почему я должен бросить знакомый C#, изучать PARI?

Вы знаете, тут больше всё-таки от задачи зависит, которую вы решаете. Ну вот, к примеру, надо вам вычислить детерминанты нескольких матриц матриц Гильберта очень большого размера или там факториал или ещё что, причём где-то параллельно. Соответственно (если не писать всё руками в С# с нуля) придётся воспользоваться подходящей математической библиотекой. Если Pari имеет всю необходимую для решения задачки математику, то можно воспользоваться ею. И тут есть несколько вариантов. Можно пользоваться Pari просто как алгебраическим калькулятором, оставаясь в рамках предложенного скриптового языка gp, тогда да, parfor и всё остальное — наше всё, вызывать это дело надо будет из командной строки, гнать вычисления в файл, тут всё просто. Либо, если надо тесно интегрировать в C# (или другую среду), то можно написать всё на Си, используя libpari.dll (кстати, конвертера gp2c я под Windows не вижу, но есть исходники), ну и подключать как динамическую библиотеку. Можно враппер написать, вот, к примеру для питона я вижу cypari/cypari2. В этом случае нам parfor вообще никаким боком не нужен (да мы его из Си/C#/Python и не можем использовать, разве что через gp_embedded_init()), тем не менее для параллелизации вычислений есть несколько возможностей (судя по примерам).

Можно воспользоваться тем, что предоставляет Pari, через mt_queue*, тут запускается три потока для параллельной обработки N1 N2 и М:

#include <pari/pari.h>

GEN
Cworker(GEN d, long kind) { return kind? det(d): Z_factor(d); }

int
main(void)
{
  long i, taskid, pending;
  GEN M,N1,N2, in,out, done;
  struct pari_mt pt;
  entree ep = {"_worker",0,(void*)Cworker,20,"GL",""};
  /* initialize PARI, postponing parallelism initialization */
  pari_init_opts(8000000,500000, INIT_JMPm|INIT_SIGm|INIT_DFTm|INIT_noIMTm);
  pari_add_function(&ep); /* add Cworker function to gp */
  pari_mt_init(); /* ... THEN initialize parallelism */
  /* Create inputs and room for output in main PARI stack */
  N1 = addis(int2n(256), 1); /* 2^256 + 1 */
  N2 = subis(int2n(193), 1); /* 2^193 - 1 */
  M = mathilbert(80);
  in  = mkvec3(mkvec2(N1,gen_0), mkvec2(N2,gen_0), mkvec2(M,gen_1));
  out = cgetg(4,t_VEC);
  /* Initialize parallel evaluation of Cworker */
  mt_queue_start(&pt, strtofunction("_worker"));
  for (i = 1; i <= 3 || pending; i++)
  { /* submit job (in) and get result (out) */
    mt_queue_submit(&pt, i, i<=3? gel(in,i): NULL);
    done = mt_queue_get(&pt, &taskid, &pending);
    if (done) gel(out,taskid) = done;
  }
  mt_queue_end(&pt); /* end parallelism */
  output(out); pari_close(); return 0;
}

Либо можно воспользоваться старым добрым openmp, наброшенным поверх pari_thread_*(), тогда это вот так будет выглядеть:

#include <pari/pari.h> /* Include PARI headers */

#include <omp.h>       /* Include OpenMP headers */

#define MAXTHREADS 3  /* Max number of parallel threads */

int
main(void)
{
  GEN M,N1,N2, F1,F2,D;
  struct pari_thread pth[MAXTHREADS];
  int numth = omp_get_max_threads(), i;
  /* Initialise the main PARI stack and global objects (gen_0, etc.) */
  pari_init(8000000,500000);
  if (numth > MAXTHREADS) {
    numth = MAXTHREADS;
    omp_set_num_threads(numth);
  }
  /* Compute in the main PARI stack */
  N1 = addis(int2n(256), 1); /* 2^256 + 1 */
  N2 = subis(int2n(193), 1); /* 2^193 - 1 */
  M = mathilbert(80);
  /*Allocate pari thread structures */
  for (i = 1; i < numth; i++) pari_thread_alloc(&pth[i],8000000,NULL);
#pragma omp parallel
  {
    int this_th = omp_get_thread_num();
    if (this_th) (void)pari_thread_start(&pth[this_th]);
#pragma omp sections
    {
#pragma omp section
      {
        F1 = factor(N1);
      }
#pragma omp section
      {
        F2 = factor(N2);
      }
#pragma omp section
      {
        D = det(M);
      }
    } /* omp sections */
    if (this_th) pari_thread_close();
  } /* omp parallel */
  pari_printf("F1=%Ps\nF2=%Ps\nlog(D)=%Ps\n", F1, F2, glog(D,3));
  for (i = 1; i < numth; i++) pari_thread_free(&pth[i]);
  return 0;
}

Ну и для любителей хардкора можно на классических POSIX потоках (эту простыню я под спойлер засуну):

Скрытый текст
#include <pari/pari.h> /* Include PARI headers */

#include <pthread.h>   /* Include POSIX threads headers */

void *
mydet(void *arg)
{
  GEN F, M;
  /* Set up thread stack and get thread parameter */
  M = pari_thread_start((struct pari_thread*) arg);
  F = QM_det(M);
  /* Free memory used by the thread */
  pari_thread_close();
  return (void*)F;
}

void *
myfactor(void *arg)  /* same principle */
{
  GEN F, N;
  N = pari_thread_start((struct pari_thread*) arg);
  F = factor(N);
  pari_thread_close();
  return (void*)F;
}

int
main(void)
{
  long prec = DEFAULTPREC;
  GEN M1,M2, N1,N2, F1,F2, D1,D2;
  pthread_t th1, th2, th3, th4; /* POSIX-thread variables */
  struct pari_thread pth1, pth2, pth3, pth4; /* pari thread variables */

  /* Initialise the main PARI stack and global objects (gen_0, etc.) */
  pari_init(32000000,500000);
  /* Compute in the main PARI stack */
  N1 = addis(int2n(256), 1); /* 2^256 + 1 */
  N2 = subis(int2n(193), 1); /* 2^193 - 1 */
  M1 = mathilbert(149);
  M2 = mathilbert(150);
  /* Allocate pari thread structures */
  pari_thread_alloc(&pth1,8000000,N1);
  pari_thread_alloc(&pth2,8000000,N2);
  pari_thread_alloc(&pth3,32000000,M1);
  pari_thread_alloc(&pth4,32000000,M2);
  /* pthread_create() and pthread_join() are standard POSIX-thread
   * functions to start and get the result of threads. */
  pthread_create(&th1,NULL, &myfactor, (void*)&pth1);
  pthread_create(&th2,NULL, &myfactor, (void*)&pth2);
  pthread_create(&th3,NULL, &mydet,    (void*)&pth3);
  pthread_create(&th4,NULL, &mydet,    (void*)&pth4); /* Start 4 threads */
  pthread_join(th1,(void*)&F1);
  pthread_join(th2,(void*)&F2);
  pthread_join(th3,(void*)&D1);
  pthread_join(th4,(void*)&D2); /* Wait for termination, get the results */
  pari_printf("F1=%Ps\nF2=%Ps\nlog(D1)=%Ps\nlog(D2)=%Ps\n",
              F1,F2, glog(D1,prec),glog(D2,prec));
  pari_thread_free(&pth1);
  pari_thread_free(&pth2);
  pari_thread_free(&pth3);
  pari_thread_free(&pth4); /* clean up */
  return 0;
}

Ну то есть при использовании Pari в многопоточной среде и интеграции со своей средой многопоточность среды остаётся как есть, она заменяет все эти "parfor", то есть если вы пользуетесь С#, то и будете использовать те же привычные Parallel.For/Parallel.Foreach, ну или я если в LabVIEW, то буду использовать свои конструкции. Я не думаю, что параллелизация parfor какая-то особенная, чудес не бывает, но всегда можно запустить бенчмарк на одинаковых данных и проверить. Как-то так.

Вы знаете, тут больше всё-таки от задачи зависит, которую вы решаете

нет, вопрос был к оформлению статьи и к вопросам, которые будут задаваться неизбежно, из статьи совсем не следует, что PARI именно такой как описываете вы. Если бы следовало, то вопросы бы не возникли.

простите, остальное читать не стал. Мы с автором эту тему закрыли, все всем довольны.

Спасибо большое за примеры параллельного кода!

Когда я компилировал с помощью gp2c, то увидел, что libpari имеет свои средства для параллелизации. А точнее итераторы для параллельного исполнения parfor_t. И функции parfor_init, parfor_stop. Но у меня не вышло их использовать без доработки.

Что нагенерил gp2c
// То, что должно работать параллельно
GEN BruteForcePart(GEN pub, GEN cnt, GEN nMaxLevel)
{
  pari_sp ltop = avma;
  setrand(gpextern("od -vAn -N8 -tu8 < /dev/urandom"));
  /* Моя секретная функция */
  avma = ltop;
  return gen_0;
}

// Для parfor_init нужно замыкание. Оно создаётся в виде функции
GEN anon_0(GEN i, GEN pub, long prec)
{
  pari_sp ltop = avma;
  GEN p1 = gen_0;
  p1 = BruteForcePart(pub, nMaxCalcInThread, nMaxX, prec);
  p1 = gerepilecopy(ltop, p1);
  return p1;
}

// В этой функции мы запускаем parfor в GP-скрипте
void transpose(GEN pub, GEN a)         /* void */
{
  pari_sp ltop = avma;
  GEN nMaxCount = gen_0;
  
  nMaxCount = gmulgs(default0("nbthreads", NULL), 100);
  {
    pari_sp btop = avma;
    while (1)
    {
      {
        parfor_t iter;    /* parfor */
        parfor_init(&iter, gen_1, nMaxCount, strtoclosure("anon_0", 1, pub));
        {
          pari_sp btop = avma;
          GEN i = gen_0, r = gen_0, p1 = gen_0;
          while ((p1 = parfor_next(&iter)))
          {
            r = gcopy(gel(p1, 2));
            i = gcopy(gel(p1, 1));
            if (gequal0(r))
              continue;
            printf0("%s %s %d %d %d %d %d %s\n", mkvecn(8, gel(a, 1), gel(a, 2), gel(r, 1), lift(gel(gel(r, 2), 1)), lift(gel(gel(r, 2), 2)), lift(gel(pub, 1)), lift(gel(pub, 2)), gel(a, 5)));
            parfor_stop(&iter);
            goto label1;
            if (gc_needed(btop, 1))
              gerepileall(btop, 3, &i, &r, &p1);
          }
        }
      }
      avma = btop;
    }
    label1:;
  }
  avma = ltop;
  return;
}

Этот код у меня не заработал. Подозреваю, что дело в том, что замыкание должно быть создано с помощью динамической компиляции GP-кода. Я не разобрался, как можно скомпилированную Си-функцию подсунуть в parfor_init.

О, спасибо огромное, а то я собрался уже было WSL расчехлять. Теперь понятнее.

Я, кстати, попробовал навскидку gp_embedded_init()/gp_embedded() и эти функции работают, что даёт возможность напрямую вызывать скрипт вообще откуда угодно. Это не самый эфффективный способ (грубо говоря там интерпретатор работает), но удобно. Офигенная игрушка, на самом деле, спасибо за наводку.

Sign up to leave a comment.