Hello, The inconsistency is still there. I pulled petsc-dev today and came across with the same problem: src/vec/vec/impls/seq/bvec2.c 239: *z = PetscRealPart(BLASdot_(&bn,xx,&one,xx,&one)); 240: *z = PetscSqrtReal(*z); This doesn't work with Intel MKL. However in one of the mailing list discussions it was said there is now check for that. This check didn't work for me apparently. Any ideas why? Regards, Alexander On 19.01.2012 17:55, Alexander Grayver wrote:
This (dirty) patch solves the problem:
diff -r a3e9ca59ab58 src/vec/vec/impls/seq/bvec2.c --- a/src/vec/vec/impls/seq/bvec2.c Tue Jan 17 22:04:05 2012 -0600 +++ b/src/vec/vec/impls/seq/bvec2.c Thu Jan 19 17:28:39 2012 +0100 @@ -232,12 +232,13 @@ PetscErrorCode ierr; PetscInt n = xin->map->n; PetscBLASInt one = 1, bn = PetscBLASIntCast(n); + PetscScalar cnorm;
PetscFunctionBegin; if (type == NORM_2 || type == NORM_FROBENIUS) { ierr = VecGetArrayRead(xin,&xx);CHKERRQ(ierr); - *z = BLASdot_(&bn,xx,&one,xx,&one); - *z = PetscSqrtReal(*z); + zdotc(&cnorm,&bn,xx,&one,xx,&one); + *z = PetscSqrtReal(PetscAbsScalar(cnorm)); ierr = VecRestoreArrayRead(xin,&xx);CHKERRQ(ierr); ierr = PetscLogFlops(PetscMax(2.0*n-1,0.0));CHKERRQ(ierr); } else if (type == NORM_INFINITY) {
The same is applied to mpi vector implementation from /petsc-dev/src/vec/vec/impls/mpi/pvec2.c Of course it works only if one uses Intel MKL BLAS/LAPACK and double complex arithmetics. Which is my case.
On 19.01.2012 15:30, Alexander Grayver wrote:
Dear petsc team,
I've been struggling with this problem for a long time and seeking for your advice now. Let's take simple petsc program:
#include <petscksp.h> int main(int argc,char **args) { Vec b; PetscReal norm;
PetscInitialize(&argc,&args,(char *)0,NULL);
VecCreateSeq(PETSC_COMM_SELF,100,&b); VecSet(b,2.0); VecNorm(b,NORM_1_AND_2,&norm);
return 0; }
This program works well if I compile petsc with non-mkl blas/lapack. However, if I compile petsc with mkl blas/lapack this program crashes:
zdotc, FP=7fff8f121e30 VecNorm_Seq, FP=7fff8f121f90 VecNorm_Seq, FP=7fff8f1220f0 VecNorm, FP=7fff8f122230 main, FP=7fff8f122290
It crashes in zdotc. You call this routine as following:
*z = BLASdot_(&bn,xx,&one,xx,&one);
When I look at the zdotc interface in mkl.h I see: void ZDOTC(MKL_Complex16 *pres, const MKL_INT *n, const MKL_Complex16 *x, const MKL_INT *incx, const MKL_Complex16 *y, const MKL_INT *incy);
I also found example here: http://software.intel.com/en-us/articles/intel-math-kernel-library-intel-mkl... There are 6 input parameters. But in "classical" BLAS implementation 5 input and 1 output.
For example, NORM_1 works well since interface to dzasum is the same.
Any ideas?
-- Regards, Alexander