Getting dgemms or dgemm to work in C

This is a discussion on Getting dgemms or dgemm to work in C within the C Programming forums, part of the General Programming Boards category; Hi there, I cannot for the life of me get dgemm to spit out the correct answer when I call ...

  1. #1
    Registered User
    Join Date
    Sep 2008
    Posts
    6

    Getting dgemms or dgemm to work in C

    Hi there,

    I cannot for the life of me get dgemm to spit out the correct answer when I call it from C.

    I would love to have a function that that was some thing like:

    void matmult(float ** A, float ** B, float ** C, int Arows, int innerdim, int Bcols)

    where you could pass it three arrays float ** or float * form (I don't care, I just want something to work) the rows of the A matrix, the shared inner dimension, and the columns of the B matrix and have the answer dumped in to the C matrix.

    Any thoughts? Figuring out how to correct for the column major order and the leading dimensions is difficult.

    Thanks so much,
    Mike

    Here's what I have so far, and it only works for an NxN times and NxN matrix:

    Code:
    /************************************
                  INCLUDES
    ************************************/
    #include <stdio.h>
    #include <stdlib.h>
    #include <math.h>
    
    /************************************
           CREATE A 1D ARRAY
    ************************************/
    float* alloc_float_1d(int ni){
      int i;
      
      float* a = (float*) malloc(ni*sizeof(float));
      
      for(i=0;i<ni;i++)
        a[i]=0;
      
      return a;
    } 
    
    
    /************************************
           1D ARRAY INDEX
    ************************************/
    int aind(int row, int col, int numcols)
    {
      return (numcols*row)+col;
    }
    
    void aprint1d(float * A, int rows, int cols)
    {
     for(int i=0;i<rows;i++)
        {
          for(int j=0;j<cols;j++)
    	printf("%1.4f ", A[aind(i,j,cols)]);
          
          printf("\n");
        }
    }
    
    
    /************************************
         Multiply 2 2d(float *) ARRAYS
    ************************************/
    void matmult(float * A, float * B, float * C, int arows, int idim, int bcols)
    {
    
      //Perform the matrix multiplication
      sgemms(A, &arows, "N", B, &bcols, "N", C, &bcols, &arows, &bcols, &idim, 0, 0);
    
    
      float temp;
      
      //FLIP UD 
      for(int i=0;i<(int)floor(arows/2.0);i++)
        {
          for(int j=0;j<bcols;j++)
    	{
    	  temp=C[aind(i,j,bcols)];
    	  C[aind(i,j,bcols)]=C[aind((bcols-1)-i,j,bcols)];
    	  C[aind((bcols-1)-i,j,bcols)]=temp;
    	}
        }
      
      //FLIP LR
      for(int i=0;i<(int)floor(arows/2.0);i++)
        {
          for(int j=0;j<bcols;j++)
    	{
    	  temp=C[aind(j,i,bcols)];
    	  C[aind(j,i,bcols)]=C[aind(j,(bcols-1)-i,bcols)];
    	  C[aind(j,(bcols-1)-i,bcols)]=temp;
    	}
    	}
    }
    
    
    //***********************************
    //***********************************
    //             MAIN
    //***********************************
    //***********************************
    void main()
    {
      //Square Array size
      int  i, j;
      int arows=3;  
      int idim=3;
      int bcols=3;
      
      //Allocate the arrays
      float * A=alloc_float_1d(arows*idim);
      float * B=alloc_float_1d(idim*bcols);
      float * C=alloc_float_1d(arows*bcols);
      
      int c=1; 
       for(i=0; i<arows; i++){
        for(j=0; j<idim; j++){ 
          A[aind(i,j,idim)] = c++;
        }
      }
      
      c=1;
      for(i=0; i<idim; i++){
        for(j=0; j<bcols; j++){ 
          B[aind(i,j,bcols)] = (idim*bcols+1)-c++;
        }
      }
      
      matmult(A,B,C,arows,idim,bcols);
       
      printf("A:\n");
      aprint1d(A,arows,idim);
      printf("\n");
      
      printf("B:\n");
      aprint1d(B,idim,bcols);
      printf("\n");
      
      printf("C:\n");
      aprint1d(C,arows,bcols);
      
    }

  2. #2
    Frequently Quite Prolix dwks's Avatar
    Join Date
    Apr 2005
    Location
    Canada
    Posts
    8,046
    void main() is a bad idea. cpwiki.sf.net/Void_Main

    What is sgemms()? You haven't defined it anywhere, nor are you including any extra header files.

    A swap(float *x, float *y) function might eliminate some duplicated code.

    [edit] Sorry I'm not of much help, my web browser is acting up . . . . [/edit]
    dwk

    Seek and ye shall find. quaere et invenies.

    "Simplicity does not precede complexity, but follows it." -- Alan Perlis
    "Testing can only prove the presence of bugs, not their absence." -- Edsger Dijkstra
    "The only real mistake is the one from which we learn nothing." -- John Powell


    Other boards: DaniWeb, TPS
    Unofficial Wiki FAQ: cpwiki.sf.net

    My website: http://dwks.theprogrammingsite.com/
    Projects: codeform, xuni, atlantis, nort, etc.

  3. #3
    Registered User
    Join Date
    Sep 2008
    Posts
    6

    dgemms

    It's a LAPACK/BLAS function that does matrix multiplication. I am using the ibm ESSL libraries. I was hoping someone had some prior experience with this particular function. I've looked on the web and it seems to give lots of people issues calling it in C since it is coded in FORTRAN which uses column-major order.

    The documentation is here:
    http://publib.boulder.ibm.com/infoce...1_hsgemms.html

    Thanks!

  4. #4
    Registered User
    Join Date
    Sep 2008
    Posts
    6
    And this is just a simple test-bed for the function. I will be using it in something much more complex.

  5. #5
    Banned master5001's Avatar
    Join Date
    Aug 2001
    Location
    Visalia, CA, USA
    Posts
    3,685
    This may sound sort of stupid but are you 100&#37; sure about your on-paper math?

  6. #6
    Registered User
    Join Date
    Sep 2008
    Posts
    6

    Yes

    It's in matlab, which doesn't lie.

  7. #7
    Registered User
    Join Date
    Sep 2008
    Posts
    6
    I just need an optimized c matrix multiplication library. Has anyone used LAPACK++ or the like?

  8. #8
    Banned master5001's Avatar
    Join Date
    Aug 2001
    Location
    Visalia, CA, USA
    Posts
    3,685
    Matlab has known limitations about its accuracy. But your point is valid enough I suppose. I don't specifically see a problem in your code, but your description of what is wrong is not exactly easy to decipher.

  9. #9
    Registered User
    Join Date
    Sep 2008
    Posts
    6
    I'm am simply looking for someone who has used fortran BLAS/LAPACK before in C and knows how to call these fortran functions correcting for the column major order.

    Thanks!!

  10. #10
    Kernel hacker
    Join Date
    Jul 2007
    Location
    Farncombe, Surrey, England
    Posts
    15,677
    It may not help much, but AMD Core Math Library's acml.h has C declarations to call BLAS and LAPACK functions. I have never used that library, so I don't know whether they are identical or just "close to" the standard BLAS/LAPACK routines.

    --
    Mats
    Compilers can produce warnings - make the compiler programmers happy: Use them!
    Please don't PM me for help - and no, I don't do help over instant messengers.

Popular pages Recent additions subscribe to a feed

Similar Threads

  1. strcmp returning 1...
    By Axel in forum C Programming
    Replies: 12
    Last Post: 09-08-2006, 07:48 PM
  2. getline() don't want to work anymore...
    By mikahell in forum C++ Programming
    Replies: 7
    Last Post: 07-31-2006, 10:50 AM
  3. Why don't the tutorials on this site work on my computer?
    By jsrig88 in forum C++ Programming
    Replies: 3
    Last Post: 05-15-2006, 10:39 PM
  4. fopen();
    By GanglyLamb in forum C Programming
    Replies: 8
    Last Post: 11-03-2002, 11:39 AM
  5. DLL __cdecl doesnt seem to work?
    By Xei in forum C++ Programming
    Replies: 6
    Last Post: 08-21-2002, 04:36 PM

Tags for this Thread


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21