How do I write my multithreading code correctly so that it is faster than a single thread? I have the following code that I am trying to optimize. I am running it on a Linux machine with 24 cores. I thought I could use multithreading to make it faster, but it's somehow making it way slower than the single thread. I believe it is because of cache misses but I am not sure how to structure this code to avoid them.

#include &lt;stdio.h&gt;              
#include &lt;stdlib.h&gt;             

#include "mex.h"
#include "string.h"
#include "pthread.h"

/* NOT sure why multi-threading is worse than single thread */
#define NTHREADS 1

/* An argument list structure */
struct argList{
    double *X;
    double *p;
    int cycles;
    int N1;
    int N2;
};
typedef struct argList aList;


long int bino_rand(double p, double n)
{
  long int    x = 0;             
  int i;
  for (i=0; i&lt;n; i++)
    if (((double)rand())/RAND_MAX &lt; p) x++;

  return(x);
}

void *bino_process(void* arglist){
    double* X   = ((aList*)arglist)-&gt;X;
    double* p   = ((aList*)arglist)-&gt;p;
    int cycles  = ((aList*)arglist)-&gt;cycles;
    int N1      = ((aList*)arglist)-&gt;N1;
    int N2      = ((aList*)arglist)-&gt;N2;
    int i,j;

/* mexPrintf("N1=%d, N2=%d, cycles=%d\n", N1, N2, cycles); */

for(i = N1;i&lt;N2;i++){
    double tmp = X[i];
    double tmp_p = p[i];
    for(j=0;j&lt;cycles;j++){
        tmp = tmp + bino_rand(tmp_p,tmp);
    }
    X[i] = tmp;
}
}


void mexFunction( int nlhs, mxArray *plhs[],
              int nrhs, const mxArray *prhs[])
{
double* n;              /* input scalar */
double* p;
int cycles;
mwSize N;               /* size of matrix */
double *X;              /* output matrix */
int i;

pthread_t threadA[NTHREADS];
int iterA[NTHREADS];
aList arglistA[NTHREADS];

int N1[NTHREADS];
int N2[NTHREADS];


/* get the value of the scalar input  */
n = mxGetPr(prhs[0]);
p = mxGetPr(prhs[1]);
cycles = (int)mxGetScalar(prhs[2]);

N = mxGetM(prhs[0]);

/* create the output matrix */
plhs[0] = mxCreateDoubleMatrix(N,1,mxREAL);

/* get a pointer to the real data in the output matrix */
X = mxGetPr(plhs[0]);

memcpy(X,n,N*sizeof(double));


for(i=0;i&lt;NTHREADS;i++){
    N1[i] = i*(N/NTHREADS);
    N2[i] = (i+1)*(N/NTHREADS);
    arglistA[i].X = X;
    arglistA[i].p = p;
    arglistA[i].cycles = cycles;


    arglistA[i].X = (double*)mxCalloc(N2[i]-N1[i],sizeof(double));
    memcpy(arglistA[i].X, &amp;X[N1[i]], (N2[i]-N1[i])*sizeof(double));

    arglistA[i].N2 = N2[i]-N1[i];        
    arglistA[i].N1 = 0;

}

/* Spawn the threads */
for(i=0;i&lt;NTHREADS;i++){
    iterA[i] = pthread_create(&amp;threadA[i],NULL,bino_process,(void*)&amp;arglistA[i]);       
}

/* Join back the threads */
for(i=0;i&lt;NTHREADS;i++){
    pthread_join(threadA[i],NULL);      
}

/* Copy back the data from the chunks into X */
for(i=0;i&lt;NTHREADS;i++){
    memcpy(&amp;X[N1[i]], arglistA[i].X, (N2[i]-N1[i])*sizeof(double));
}

for(i=0;i&lt;NTHREADS;i++){
    mxFree(arglistA[i].X);
}

}
How did you time your program?
I ran a modified version of your program (I didn't have &quot;mex.h&quot; so I rewrote the parts that needed it).  When I ran your program with 1, 2, and 4 threads, I got the following results:
1 thread  = 9.38 sec
2 threads = 4.93 sec
4 threads = 2.81 sec

So it seems to be working for me.  How did you time your program?  You may have been using a function such as clock() that measures cpu time across all cores, instead of measuring wall clock time.  See this Stackoverflow question for more information on that.