Я сначала переписал его, чтобы было немного легче следовать:
#include <math.h>
double Param_RE_Tterm_approx(double Tterm, double const* par) {
double value = 0.;
if (Tterm > 0.) {
// time after Che angle peak
if ( fabs(Tterm/ par[0]) >= 1.e-2) {
value += -1./(par[0])*exp(-1.*Tterm/(par[0]));
} else {
value += -1./par[0]*(1. - Tterm/par[0] + Tterm*Tterm/(par[0]*par[0]*2.) - Tterm*Tterm*Tterm/(par[0]*par[0]*par[0]*6.) );
}
if ( fabs(Tterm* par[1]) >= 1.e-2) {
value += par[2]* par[1]*pow( 1.+par[1]*Tterm, par[2]-1. );
} else {
value += par[2]*par[1]*( 1.+(par[2]-1.)*par[1]*Tterm + (par[2]-1.)*(par[2]-1.-1.)/2.*par[1]*par[1]*Tterm*Tterm + (par[2]-1.)*(par[2]-1.-1.)*(par[2]-1.-2.)/6.*par[1]*par[1]*par[1]*Tterm*Tterm*Tterm );
}
} else {
// time before Che angle peak
if ( fabs(Tterm/ par[3]) >= 1.e-2 ) {
value += -1./ par[3] *exp(-1.*Tterm/ par[3]);
} else {
value += -1./par[3]*(1. - Tterm/par[3] + Tterm*Tterm/(par[3]*par[3]*2.) - Tterm*Tterm*Tterm/(par[3]*par[3]*par[3]*6.) );
}
if ( fabs(Tterm* par[4]) >= 1.e-2 ) {
value += par[5]* par[4] *pow( 1.+ par[4]*Tterm, par[5]-1. );
} else {
value += par[5]*par[4]*( 1.+(par[5]-1.)*par[4]*Tterm + (par[5]-1.)*(par[5]-1.-1.)/2.*par[4]*par[4]*Tterm*Tterm + (par[5]-1.)*(par[5]-1.-1.)*(par[5]-1.-2.)/6.*par[4]*par[4]*par[4]*Tterm*Tterm*Tterm );
}
}
return value * 1.e9;
}
Затем мы можем взглянуть на его структуру.
Существует две основные ветви - Tterm отрицательная (до) и положительная (после). Они соответствуют использованию 0,1,2 или 3,4,5 в массиве par.
Тогда в каждом случае мы делаем две вещи, чтобы увеличить стоимость. В обоих случаях для небольших случаев мы используем полином, а для больших случаев - уравнение экспоненты / степени.
Как и догадка , это связано с тем, что многочлен является приличным приближением для экспоненты для малых значений - ошибка приемлема. Что вы должны сделать, это подтвердить эту догадку - взгляните на разложение в ряд Тейлора "большого" уравнения, основанного на степени / экспоненте, и посмотрите, как оно согласуется с полиномами. Или проверить численно.
Если это так, это означает, что в этом уравнении есть известная величина ошибки, которая является приемлемой. Довольно часто существуют более быстрые версии exp или pow, которые имеют известную величину максимальной ошибки; подумайте об их использовании.
Если это не так, все равно может быть допустимое количество ошибок, но приближение ряда Тейлора может дать вам информацию «в коде» о том, что является приемлемым количеством ошибок.
Следующий шаг, который я предприму, - это разорвать 8 частей этого уравнения на части. Существует положительный / отрицательный, первый и второй value+= в каждой ветви, а затем полиномиальный / экспоненциальный случай.
Я догадываюсь, что exp занимает ~ 1/3 времени pow, потому что у вас есть 3 вызова pow к 1 вызову exp в вашей функции, но вы можете найти что-то интересное типа "все наши время на самом деле в Tterm> 0. case "или что у вас.
Теперь рассмотрим сайты вызовов. Есть ли в Tterm шаблон, по которому вы передаете эту функцию? Т.е. вы склонны сдавать Tterms в грубо отсортированном порядке? Если это так, вы можете выполнить проверку, для какой функции вызывать за пределами вызова этой функции, и делать это партиями.
Простое выполнение этого в пакетном режиме и компиляция с оптимизацией и вставкой тел функций может привести к удивительной разнице; компиляторы становятся лучше при векторизации работы.
Если это не сработает, вы можете начать что-то делать. На современном компьютере вы можете иметь 4-60 потоков, решающих эту проблему независимо, и эта проблема выглядит так, как будто бы вы получили почти линейное ускорение. Для такой задачи подойдет базовая библиотека потоков, такая как TBB.
Для следующего шага, если вы получаете большие партии данных и вам нужно много обработать, вы можете поместить их в графический процессор и там решить. К сожалению, связь с GPU <-> RAM невелика, поэтому простое выполнение математических операций в этой функции на GPU и чтение / запись назад и вперед с RAM не даст вам много, если вообще произойдет. Но если на GPU может работать больше, чем просто это, это может стоить.