daxpy dcopy计算

daxpy每秒计算次数:

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include "cblas.h"

#define ARRAY_LENGTH 1000000

void main() {


    int n;                          /*! array size */
    double da;                      /*! double constant */
    double *dx;                     /*! input double array */
    int incx;                        /*! input stride */
    double *dy;                      /*! output double array */
    int incy;                       /*! output stride */

    int i;

    int warp_count = 0;
    int max_warp = 1000;
    long int count = 0;
    time_t b_second,l_second;

    time_t rawtime;
    struct tm * timeinfo;

    n = ARRAY_LENGTH;
    da = 99999;
    dx = (double*)malloc(sizeof(double)*n);
    incx = 1;
    dy = (double*)malloc(sizeof(double)*n);
    incy = 1;

    srand((unsigned)time(NULL));

    for(i=0;i<n;i++){
        dx[i] = rand();
        dy[i] = rand();
        //printf("%f ",dy[i]);    //输出原来的dy
    }

    while(1){

        b_second = time(NULL);
        l_second = b_second+1;

        while ((b_second=time(NULL))<l_second) {
            cblas_daxpy(n, da, dx,incx, dy, incy);  //运行daxpy程序
            count++;
        }

        time(&rawtime);
        timeinfo = localtime (&rawtime);
        printf("Time: %s ", asctime (timeinfo));
        printf("%ld\n",count);

        count=0;
        warp_count++;
        if(warp_count==max_warp){
            break;
        }
    }
}

dcopy每秒计算次数:

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include "cblas.h"

#define ARRAY_LENGTH 100000

void main() {


    int n;                          /*! array size */
    double *dx;                     /*! input double array */
    int incx;                        /*! input stride */
    double *dy;                      /*! output double array */
    int incy;                       /*! output stride */

    int i;

    int warp_count = 0;
    int max_warp = 1000;
    long int count = 0;
    time_t b_second,l_second;

    time_t rawtime;
    struct tm * timeinfo;

    n = ARRAY_LENGTH;
    dx = (double*)malloc(sizeof(double)*n);
    incx = 1;
    dy = (double*)malloc(sizeof(double)*n);
    incy = 1;

    srand((unsigned)time(NULL));

    for(i=0;i<n;i++){
        dx[i] = rand();
        dy[i] = rand();
        //printf("%f ",dy[i]);    //输出原来的dy
    }

    while(1){

        b_second = time(NULL);
        l_second = b_second+1;

        while ((b_second=time(NULL))<l_second) {
            cblas_dcopy(n, dx,incx, dy, incy);      //运行dcopy程序
            count++;
        }

        time(&rawtime);
        timeinfo = localtime (&rawtime);
        printf("Time: %s ", asctime (timeinfo));
        printf("%ld\n",count);

        count=0;
        warp_count++;
        if(warp_count==max_warp){
            break;
        }
    }
}
subroutine FGMRES(A, IA, JA, RHS, N, nzmax, COMPUTED_SOLUTION) implicit none integer N integer SIZE parameter (SIZE=128) integer, intent(in) :: nzmax !--------------------------------------------------------------------------- ! Define arrays for the upper triangle of the coefficient matrix ! Compressed sparse row storage is used for sparse representation !--------------------------------------------------------------------------- integer IA(N+1) integer JA(nzmax) real*8 A(nzmax) !--------------------------------------------------------------------------- ! Allocate storage for the ?par parameters and the solution/rhs/residual vectors !--------------------------------------------------------------------------- integer IPAR(SIZE) real*8 DPAR(SIZE), TMP(N*(2*150+1)+(150*(150+9))/2+1) real*8 RHS(N) !, B(N) real*8 COMPUTED_SOLUTION(N) ! real*8 RESIDUAL(N) !--------------------------------------------------------------------------- ! Some additional variables to use with the RCI (P)FGMRES solver !--------------------------------------------------------------------------- integer ITERCOUNT, EXPECTED_ITERCOUNT parameter (EXPECTED_ITERCOUNT=233) integer RCI_REQUEST, I real*8 DVAR !--------------------------------------------------------------------------- ! An external BLAS function is taken from MKL BLAS to use ! with the RCI (P)FGMRES solver ! Save the right-hand side in vector B for future use !--------------------------------------------------------------------------- ! call DCOPY(N, RHS, 1, B, 1) !--------------------------------------------------------------------------- ! Initialize the initial guess !--------------------------------------------------------------------------- do I=1,N COMPUTED_SOLUTION(I)=0. END DO !--------------------------------------------------------------------------- ! Initialize the solver !--------------------------------------------------------------------------- CALL DFGMRES_INIT(N, COMPUTED_SOLUTION, RHS, RCI_REQUEST, IPAR, DPAR, TMP) IF (RCI_REQUEST.NE.0) GOTO 999 !--------------------------------------------------------------------------- ! Set the desired parameters: ! do the restart after 2 iterations: IPAR(15)=2 ! LOGICAL parameters: ! do not do the stopping test for the maximal number of iterations: IPAR(8)=0 ! do the Preconditioned iterations of FGMRES method: IPAR(11)=1 ! DOUBLE PRECISION parameters ! set the relative tolerance to 1.0D-3 instead of default value 1.0D-6: DPAR(1)=1.0D-3 !--------------------------------------------------------------------------- IPAR(9)=1 IPAR(10)=0 IPAR(12)=1 DPAR(1)=1.0D-6 !--------------------------------------------------------------------------- ! Check the correctness and consistency of the newly set parameters !--------------------------------------------------------------------------- CALL DFGMRES_CHECK(N, COMPUTED_SOLUTION, RHS, RCI_REQUEST, IPAR, DPAR, TMP) IF (RCI_REQUEST.NE.0) GOTO 999 !--------------------------------------------------------------------------- ! Compute the solution by RCI (P)FGMRES solver with preconditioning ! Reverse Communication starts here !--------------------------------------------------------------------------- 1 CALL DFGMRES(N, COMPUTED_SOLUTION, RHS, RCI_REQUEST, IPAR, DPAR, TMP) !--------------------------------------------------------------------------- ! If RCI_REQUEST=0, then the solution was found with the required precision !--------------------------------------------------------------------------- IF (RCI_REQUEST.EQ.0) GOTO 3 !--------------------------------------------------------------------------- ! If RCI_REQUEST=1, then compute the vector A*TMP(IPAR(22)) ! and put the result in vector TMP(IPAR(23)) !--------------------------------------------------------------------------- IF (RCI_REQUEST.EQ.1) THEN CALL MKL_DCSRGEMV('N',N, A, IA, JA, TMP(IPAR(22)), TMP(IPAR(23))) GOTO 1 ! END IF !!--------------------------------------------------------------------------- !! If RCI_request=2, then do the user-defined stopping test !! The residual stopping test for the computed solution is performed here !!--------------------------------------------------------------------------- !! NOTE: from this point vector B(N) is no longer containing the right-hand !! side of the problem! It contains the current FGMRES approximation to the !! solution. If you need to keep the right-hand side, save it in some other !! vector before the call to DFGMRES routine. Here we saved it in vector !! RHS(N). The vector B is used instead of RHS to preserve the original !! right-hand side of the problem and guarantee the proper restart of FGMRES !! method. Vector B will be altered when computing the residual stopping !! criterion! !!--------------------------------------------------------------------------- ! IF (RCI_REQUEST.EQ.2) THEN !! Request to the DFGMRES_GET routine to put the solution into B(N) via IPAR(13) ! IPAR(13)=1 !! Get the current FGMRES solution in the vector B(N) ! CALL DFGMRES_GET(N, COMPUTED_SOLUTION, B, RCI_REQUEST, IPAR, DPAR, TMP, ITERCOUNT) !! Compute the current true residual via MKL (Sparse) BLAS routines ! CALL MKL_DCSRGEMV('N', N, A, IA, JA, B, RESIDUAL) ! CALL DAXPY(N, -1.0D0, RHS, 1, RESIDUAL, 1) ! DVAR=DNRM2(N, RESIDUAL, 1) ! IF (DVAR.LT.1.0E-3) THEN ! GOTO 3 ! ELSE ! GOTO 1 ! END IF !!--------------------------------------------------------------------------- !! If RCI_REQUEST=anything else, then DFGMRES subroutine failed !! to compute the solution vector: COMPUTED_SOLUTION(N) !!--------------------------------------------------------------------------- ELSE ! call data WRITE(*,*) 'This example FAILED as the solver has returned the ERROR code',RCI_REQUEST pause GOTO 999 END IF !--------------------------------------------------------------------------- ! Reverse Communication ends here ! Get the current iteration number and the FGMRES solution (DO NOT FORGET to ! call DFGMRES_GET routine as COMPUTED_SOLUTION is still containing ! the initial guess!) !--------------------------------------------------------------------------- 3 CALL DFGMRES_GET(N, COMPUTED_SOLUTION, RHS, RCI_REQUEST, IPAR, DPAR, TMP, ITERCOUNT) IF (ITERCOUNT<=EXPECTED_ITERCOUNT) THEN ! WRITE(*,*) 'This example has successfully PASSED through all steps of computation!' ! print*,'number of iteration',ITERCOUNT RETURN ELSE WRITE(*,*) 'This example have FAILED as the number of iterations differs from the expected number of iterations' pause STOP 15 END IF 999 WRITE(*,*) 'This example FAILED as the solver has returned the ERROR code',RCI_REQUEST CALL MKL_FREE_BUFFERS STOP 16 END SUBROUTINE FGMRES详细解释下这段代码,中文注释
最新发布
07-02
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值