Making a fitfunction multithreaded

Sometimes you will want to make your fit function multithreaded. This is not much use when you use Funcfit or Curvefit, as the operation does the threading itself. However, it can be very useful when using GenCurvefit as this XOP operation cannot perform multithreading itself, so you need to do it in the fitfunction. The code list below is a demonstration of how to do this. Simply change the code in the _kernel functions and you're off. You can of course change all the function names to your hearts content.

You should note that significant speedups are only obtained when the fitfunction to be calculated takes a long time. If it is too quick then the thread overhead will kill you. Test the speedups before you use code exclusively

Function myfitfunc_wrapper(w, yy, xx) : fitfunc
Wave w, yy, xx 
    Multithread yy = myfitfunc_kernel(w, xx)
End

Threadsafe function myfitfunc_kernel(w, xx)
Wave w
variable xx
    //your code goes here
    return w[0] + sqrt(xx) * w[1] * sin(xx) * cos(x) * ln(x)
End

-OR-
Function myfitfunc_AAO_wrapper(w, yy, xx) : fitfunc
Wave w, yy, xx
 
    variable nt = ThreadProcessorCount
    variable tgID = ThreadGroupCreate(nt)
    variable ii, startP, perT, remainingP
 
    perT = floor(dimsize(yy, 0) / nt)
    remainingP = dimsize(yy, 0)
 
    make/n=(nt)/free/wave theDataY, theDataX
    for(ii = 0 ; ii < nt ; ii+=1)
        if(ii < nt -1)
            theDataY[ii] = newfreewave(0x04, perT)
            theDataX[ii] = newfreewave(0x04, perT)
        else
            theDataY[ii] = newfreewave(0x04, remainingP)
            theDataX[ii] = newfreewave(0x04, remainingP)       
        endif
        Wave yyy = theDataY[ii]
        Wave xxx = theDataX[ii]
 
        yyy = yy[p+startP]
                xxx = xx[p+startP]
 
        startP += perT
        remainingP -= perT
        ThreadStart tgID, ii, myfitfunc_AAO_kernel(w, theDataY[ii], theDataX[ii])
        WaveClear yyy, xxx
    endfor
    variable finished = ThreadGroupWait(tgID, inf)
    finished = threadgrouprelease(tgID)
   
    //now copy back into structure wave
    startP = 0
    perT = floor(dimsize(yy, 0) / nt)
    remainingP = dimsize(yy, 0)
    for(ii = 0 ; ii < nt ; ii+=1)
        Wave yyy = theDataY[ii]
        yy[startP, startP+numpnts(yyy)] = yyy[p-startP]
        Waveclear yyy

        startP += perT
        remainingP -= perT
    endfor
End

Threadsafe function myfitfunc_AAO_kernel(w, yy, xx)
Wave w, yy, xx
    //your code goes here
    yy = w[0] + w[1] * sqrt(xx) * sin(xx) * cos(x) * ln(x)
End


-OR- if you wanted to do a structure fit

Structure fitfuncStruct  
    Wave w
    wave y
    wave x[50]
 
    int16 numVarMD
    wave ffsWaves[50]
    wave ffsTextWaves[10]
    variable ffsvar[5]
    string ffsstr[5]
    nvar ffsnvars[5]
    svar ffssvars[5]
    funcref allatoncefitfunction ffsfuncrefs[10]
    uint32 ffsversion    // Structure version.
    EndStructure

Threadsafe Function allatoncefitfunction(w,y,x)
    Wave w,y,x
    //you need this to act as a signature for the fitfuncStruct
    //it needs to be threadsafe
End

Threadsafe function myfitfunc_AAO_kernel(w, yy, xx)
Wave w, yy, xx
    //your code goes here
    yy = w[0] + w[1] * sqrt(xx) * sin(xx) * cos(x) * ln(x)
End

Function myfitfunc_AAO_wrapper(s) : fitfunc
    Struct fitfuncStruct &s
 
    //work out how many points per thread
    variable nt =  ThreadProcessorCount
    variable tgID = ThreadGroupCreate(nt)
    variable ii, startP, perT, remainingP, finished
    perT = floor(dimsize(s.y, 0) / nt)
    remainingP = dimsize(s.y, 0)
       
    //create waveref waves to hold references to the y and x waves to be distributed amongst each thread.
    //you could create more than y and x if your function took more arguments
    make/n=(nt)/free/wave theDataY, theDataX
   
    //create the coefficients wave to be sent to all the threads, as you can't pass in s.w to a thread
    Wave theCoefs = s.w

    //this is the function we are going to pass the coefficients, and y and x waves onto.
     //-It doesn't have to be an allatoncefitfunction it can be any function type you want (so long as it doesn't take structures).
    //- It needs to be threadsafe!
    //- the function can take as many parameters as you want, but you'll have to create a free wave for each of them (like theCoefs) to
    //    be passed into the thread.
    Funcref allatoncefitfunction theFitFunction = s.ffsfuncrefs[0]
   
    for(ii = 0 ; ii < nt ; ii+=1)
        if(ii < nt -1)
            theDataY[ii] = newfreewave(0x04, perT)
            theDataX[ii] = newfreewave(0x04, perT)
        else
            theDataY[ii] = newfreewave(0x04, remainingP)
            theDataX[ii] = newfreewave(0x04, remainingP)       
        endif
        Wave yyy = theDataY[ii]
        Wave xxx = theDataX[ii]
       
        yyy = s.y[p+startP]
        xxx = s.x[0][p+startP]
       
        startP += perT
        remainingP -= perT

        //NOTE that you can't pass in a structure parameter in a threadstart operation. See:
        //DisplayHelptopic "Threadstart"
        // for more information
        ThreadStart tgID, ii, myfitfunc_AAO_kernel(theCoefs, yyy, xxx)
        WaveClear yyy, xxx
    endfor
     finished = ThreadGroupWait(tgID, inf)
    finished = threadgrouprelease(tgID)

    //now copy back into structure wave
    startP = 0
    perT = floor(dimsize(s.y, 0) / nt)
    remainingP = dimsize(s.y, 0)
    for(ii = 0 ; ii < nt ; ii+=1)
        Wave yyy = theDataY[ii]
        s.y[startP, startP+numpnts(yyy)] = yyy[p-startP]
        Waveclear yyy

        startP += perT
        remainingP -= perT
    endfor
End
Hi, first of all thanks for your post! I am using pretty massive all-at-once fitting functions, and I wanted to see how they perform using genetic optimisation code. I have an i7 CPU, so I want to exploit multithreading. I tried your code but it didn't work, as I received the error "FitFunction returned NaN/INF for one or more values". I specify that these fitfunctions in their "normal" form are working flawlessly. Do you know the possible source of error? I attach a .txt of my code. Let me know if you need further info.
Thanks!

Riccardo
fitfunc_code.txt
can you post an experiment with some data?

Can't promise I can solve your problem, but I'll have a look.
I'll shrink the one I am using right now, to make it lighter. As soon as I can I will post it!

Forum

Support

Gallery

Igor Pro 8

Learn More

Igor XOP Toolkit

Learn More

Igor NIDAQ Tools MX

Learn More