[petsc-dev] Fwd: [petsc-maint #56111] petsc4py does not work with complex numbers?
Satish Balay
balay at mcs.anl.gov
Tue Nov 16 12:41:42 CST 2010
Barry,
I'm getting the following in the debugger:
Breakpoint 1, KSPSolve_BCGSL (ksp=0x2cde280) at bcgsl.c:144
144 if (ksp->reason) {
(gdb) p ksp->reason
$2 = KSP_DIVERGED_DTOL
And the following change fixes it for me. Does this look ok?
Satish
------------
asterix:/home/balay/tmp/petsc-dev>hg diff
diff -r 9773839bae8f src/ksp/ksp/impls/bcgsl/bcgsl.c
--- a/src/ksp/ksp/impls/bcgsl/bcgsl.c Tue Nov 16 08:39:10 2010 -0600
+++ b/src/ksp/ksp/impls/bcgsl/bcgsl.c Tue Nov 16 12:41:10 2010 -0600
@@ -146,6 +146,7 @@
ksp->its = k+j;
ksp->rnorm = nrm0;
ierr = PetscObjectGrantAccess(ksp);CHKERRQ(ierr);
+ if (ksp->reason <0) bBombed = PETSC_TRUE;
break;
}
}
asterix:/home/balay/tmp/petsc-dev>
On Tue, 16 Nov 2010, Barry Smith wrote:
>
> Satish,
>
> Check the values in the part of the code:
>
> /* Symmetrize MZa */
> for(i = 0; i <= bcgsl->ell; ++i) {
> for(j = i+1; j <= bcgsl->ell; ++j) {
> MZa[i*ldMZ+j] = MZa[j*ldMZ+i] = PetscConj(MZa[j*ldMZ+i]);
> }
> }
>
> Are all the inputs and outputs to this chunk of code reasonable?
>
> Thanks
>
> Barry
>
>
> On Nov 16, 2010, at 8:23 AM, Satish Balay wrote:
>
> > Yes - valgrind doesn't show anything with this code.. [it shows a
> > bunch of python stuff - before PetscInitialize() - so I'm guessing
> > this is unreleated python stuff..]
> >
> > And suprisingly fp_trap isn't catching the SIGFPE - when run inside
> > valgrind..
> >
> > BTW: due to the SIGFPE - there is a Nan in the var - so this check
> > [-(-var)==var] and it fails and gives the following message -
> > eventhough its sequential run.
> >
> >>>> [0] Scalar value must be same on all processes, argument # 3
> >
> > Satish
> >
> >
> > On Tue, 16 Nov 2010, Barry Smith wrote:
> >
> >>
> >> Satish,
> >>
> >> Valgrind ok?
> >>
> >> Barry
> >>
> >> On Nov 16, 2010, at 7:23 AM, Satish Balay wrote:
> >>
> >>> With '-fp_trap -start_in_debugger' in ~/.petscrc, I get a crash at:
> >>>
> >>>>>>>>>>>>
> >>> Program received signal SIGFPE, Arithmetic exception.
> >>> 0x00007f86e5ee17b4 in ATL_zreftrsmLLNN () from /usr/lib64/atlas/libatlas.so.3
> >>> (gdb) where
> >>> #0 0x00007f86e5ee17b4 in ATL_zreftrsmLLNN ()
> >>> from /usr/lib64/atlas/libatlas.so.3
> >>> #1 0x00007f86e5efc9fb in ATL_ztrsmLLNN () from /usr/lib64/atlas/libatlas.so.3
> >>> #2 0x00007f86e5efc4e1 in ATL_ztrsm () from /usr/lib64/atlas/libatlas.so.3
> >>> #3 0x000000394422ebbc in ATL_zpotrs () from /usr/lib64/atlas/liblapack.so.3
> >>> #4 0x00000039446a69cb in zpotrs_ () from /usr/lib64/atlas/liblapack.so.3
> >>> #5 0x00007f86e71df480 in KSPSolve_BCGSL (ksp=0x1372da0) at bcgsl.c:179
> >>> #6 0x00007f86e74dbd1e in KSPSolve (ksp=0x1372da0, b=0x13a9c30, x=0x13a5b20)
> >>> at itfunc.c:427
> >>> #7 0x00007f86e7ef5898 in __pyx_pf_8petsc4py_5PETSc_3KSP_solve (
> >>> __pyx_v_self=0xc8e090, __pyx_args=<value optimized out>,
> >>> __pyx_kwds=<value optimized out>) at src/petsc4py.PETSc.c:89419
> >>>
> >>>
> >>> #5 0x00007f86e71df480 in KSPSolve_BCGSL (ksp=0x1372da0) at bcgsl.c:179
> >>> 179 LAPACKpotrs_("Lower", &bell, &ione, &MZa[1+ldMZ], &ldMZ, &AY0c[1], &ldMZ, &bierr);
> >>> (gdb) p ldMZ
> >>> $1 = 3
> >>> (gdb) p ldMZ
> >>> $2 = 3
> >>> (gdb) p MZa[1+ldMZ]
> >>> $3 = 2.1147927874194505e-14 + 0 * I
> >>> (gdb) p MZa[1+ldMZ+1]
> >>> $4 = 2.1147927874194505e-14 + 0 * I
> >>> (gdb) p MZa[1+ldMZ+2]
> >>> $5 = 4.472348533721329e-28 + -0 * I
> >>> (gdb) p AY0c[1]
> >>> $6 = 2.1147927874194505e-14 + 0 * I
> >>> (gdb) p AY0c[2]
> >>> $7 = 0 + 0 * I
> >>> (gdb) p AY0c[3]
> >>> $8 = 1.9966503766457982e-314 + 9.5700515599449456e-321 * I
> >>> (gdb)
> >>> <<<<<<<<<<<<<<<
> >>>
> >>> So - perhaps some of these values are garbage - generated in earler
> >>> part of the code? I don't know where the bug is yet..
> >>>
> >>> [BTW: VecEqual fix is pushed to both petsc-31 and petsc-dev]
> >>>
> >>> Satish
> >>>
> >>> On Mon, 15 Nov 2010, Lisandro Dalcin wrote:
> >>>
> >>>> petsc-dev has a bug in VecEqual for the case of complex scalars. I've
> >>>> mailed Satish about this with a single-line patch. After applying the
> >>>> patch, almost all petsc4py tests pass.
> >>>>
> >>>> Still, there is a failure, but I doubt it is petsc4py's fault:
> >>>>
> >>>> ======================================================================
> >>>> ERROR: testSolve (test_ksp.TestKSPBCGSL)
> >>>> ----------------------------------------------------------------------
> >>>> Traceback (most recent call last):
> >>>> File "test/test_ksp.py", line 132, in testSolve
> >>>> self.ksp.solve(b, x)
> >>>> File "KSP.pyx", line 326, in petsc4py.PETSc.KSP.solve
> >>>> (src/petsc4py.PETSc.c:89520)
> >>>> Error: error code 62
> >>>> [0] KSPSolve() line 427 in /usr/local/petsc/dev/src/ksp/ksp/interface/itfunc.c
> >>>> [0] KSPSolve_BCGSL() line 241 in
> >>>> /usr/local/petsc/dev/src/ksp/ksp/impls/bcgsl/bcgsl.c
> >>>> [0] VecMAXPY() line 1202 in /usr/local/petsc/dev/src/vec/vec/interface/rvector.c
> >>>> [0] Invalid argument
> >>>> [0] Scalar value must be same on all processes, argument # 3
> >>>>
> >>>>
> >>>> On 15 November 2010 20:05, Satish Balay <balay at mcs.anl.gov> wrote:
> >>>>> Ok. For now - I've added complex=1 for petsc4py.py [and mpi4py.py]
> >>>>>
> >>>>> Satish
> >>>>>
> >>>>> On Mon, 15 Nov 2010, Lisandro Dalcin wrote:
> >>>>>
> >>>>>> On 15 November 2010 19:01, Satish Balay <balay at mcs.anl.gov> wrote:
> >>>>>>> I get the following with complex enabled. Does this mean 'some tests
> >>>>>>> run but others fail'?
> >>>>>>>
> >>>>>>> I guess we'll have to wait for Lisandro to reply..
> >>>>>>>
> >>>>>>
> >>>>>> Complex numbers should work just fine, at least with 3.1, I'm going to
> >>>>>> check against petsc-dev.
> >>>>>>
> >>>>>>>
> >>>>>>> asterix:/home/balay/tmp/petsc-dev/externalpackages/petsc4py-dev>make test
> >>>>>>> python test/runtests.py < /dev/null
> >>>>>>> [0 at asterix] Python 2.7 (/usr/bin/python)
> >>>>>>> [0 at asterix] PETSc 3.1.0-p5 development (conf: 'asterix64')
> >>>>>>> [0 at asterix] petsc4py 1.1.2 (build/lib.linux-x86_64-2.7/petsc4py)
> >>>>>>> .......................................................................................................................................................................................................................EE...................................................................................................F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F......
..
> > F
> >>> ...
> >>>>> ....
> >>>>>>> .F........F........F........F.........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F........F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F..F...................FE..........FF.FEF...........F.....................................................................................................................................................................................................................FFFF........FFFF................................................................F...F..F.......F...F..F.......F...F..F........
> >>>>>>> ======================================================================
> >>>>>>> ERROR: testSetMonitor (test_ksp.TestKSPBCGSL)
> >>>>>>> ----------------------------------------------------------------------
> >>>>>>> Traceback (most recent call last):
> >>>>>>> File "test/test_ksp.py", line 154, in testSetMonitor
> >>>>>>> self.testSolve()
> >>>>>>> File "test/test_ksp.py", line 132, in testSolve
> >>>>>>> self.ksp.solve(b, x)
> >>>>>>> File "KSP.pyx", line 326, in petsc4py.PETSc.KSP.solve (src/petsc4py.PETSc.c:89419)
> >>>>>>> Error: error code 62
> >>>>>>> [0] KSPSolve() line 427 in src/ksp/ksp/interface/itfunc.c
> >>>>>>> [0] KSPSolve_BCGSL() line 241 in src/ksp/ksp/impls/bcgsl/bcgsl.c
> >>>>>>> [0] VecMAXPY() line 1202 in src/vec/vec/interface/rvector.c
> >>>>>>> [0] Invalid argument
> >>>>>>> [0] Scalar value must be same on all processes, argument # 3
> >>>>>>>
> >>>>>>> ======================================================================
> >>>>>>> ERROR: testSolve (test_ksp.TestKSPBCGSL)
> >>>>>>> ----------------------------------------------------------------------
> >>>>>>> Traceback (most recent call last):
> >>>>>>> File "test/test_ksp.py", line 132, in testSolve
> >>>>>>> self.ksp.solve(b, x)
> >>>>>>> File "KSP.pyx", line 326, in petsc4py.PETSc.KSP.solve (src/petsc4py.PETSc.c:89419)
> >>>>>>> Error: error code 62
> >>>>>>> [0] KSPSolve() line 427 in src/ksp/ksp/interface/itfunc.c
> >>>>>>> [0] KSPSolve_BCGSL() line 241 in src/ksp/ksp/impls/bcgsl/bcgsl.c
> >>>>>>> [0] VecMAXPY() line 1202 in src/vec/vec/interface/rvector.c
> >>>>>>> [0] Invalid argument
> >>>>>>> [0] Scalar value must be same on all processes, argument # 3
> >>>>>>>
> >>>>>>> ======================================================================
> >>>>>>> ERROR: testCreateTranspose (test_mat_py.TestDiagonal)
> >>>>>>> ----------------------------------------------------------------------
> >>>>>>> Traceback (most recent call last):
> >>>>>>> File "test/test_mat_py.py", line 100, in tearDown
> >>>>>>> self.assertEqual(getrefcount(ctx), 2)
> >>>>>>> AssertionError: 3 != 2
> >>>>>>>
> >>>>>>> ======================================================================
> >>>>>>> ERROR: testMultTransposeNewMeth (test_mat_py.TestIdentity)
> >>>>>>> ---------------------------
> >>>>>>>
> >>>>>>>
> >>>>>>> On Mon, 15 Nov 2010, Barry Smith wrote:
> >>>>>>>
> >>>>>>>>
> >>>>>>>>
> >>>>>>>>>
> >>>>>>>>>
> >>>>>>>>> Could someone who knows answer this? If it supports complex then we'll need to change the flag in petsc4py.py
> >>>>>>>>>
> >>>>>>>>> Barry
> >>>>>>>>>
> >>>>>>>>> *******************************************************************************
> >>>>>>>>> UNABLE to CONFIGURE with GIVEN OPTIONS (see configure.log for details):
> >>>>>>>>> -------------------------------------------------------------------------------
> >>>>>>>>> Cannot use petsc4py with complex numbers it is not coded for this capability
> >>>>>>>>> *******************************************************************************
> >>>>>>>>> File "./config/configure.py", line 260, in petsc_configure
> >>>>>>>>> framework.configure(out = sys.stdout)
> >>>>>>>>> File "/Users/srinath/work/Packages/rev-packages/petsc-dev/config/BuildSystem/config/framework.py", line 977, in configure
> >>>>>>>>> child.configure()
> >>>>>>>>> File "/Users/srinath/work/Packages/rev-packages/petsc-dev/config/BuildSystem/config/package.py", line 488, in configure
> >>>>>>>>> self.consistencyChecks()
> >>>>>>>>> File "/Users/srinath/work/Packages/rev-packages/petsc-dev/config/PETSc/package.py", line 40, in consistencyChecks
> >>>>>>>>> raise RuntimeError('Cannot use '+self.name+' with complex numbers it is not coded for this capability')
> >>>>>>>>>
> >>>>>>>>>
> >>>>>>>>>
> >>>>>>>>> Begin forwarded message:
> >>>>>>>>>
> >>>>>>>>>> From: Srinath Vadlamani <srinath at txcorp.com>
> >>>>>>>>>> Date: November 15, 2010 3:27:55 PM CST
> >>>>>>>>>> To: petsc-maint <petsc-maint at mcs.anl.gov>
> >>>>>>>>>> Subject: [petsc-maint #56111] petsc4py does not work with complex numbers?
> >>>>>>>>>> Reply-To: petsc-maint at mcs.anl.gov, Srinath Vadlamani <srinath at txcorp.com>
> >>>>>>>>>>
> >>>>>>>>>> Is this true?
> >>>>>>>>>>
> >>>>>>>>>> --
> >>>>>>>>>> =========================
> >>>>>>>>>> Srinath Vadlamani, PhD.
> >>>>>>>>>> Tech-X Corp. Research Scientist
> >>>>>>>>>> 5621 Arapahoe Ave. Suite A
> >>>>>>>>>> Boulder, CO 80303
> >>>>>>>>>> 303-996-2034
> >>>>>>>>>> =========================
> >>>>>>>>>>
> >>>>>>>>>> Is this true?
> >>>>>>>>>>
> >>>>>>>>>> --
> >>>>>>>>>> =========================
> >>>>>>>>>> Srinath Vadlamani, PhD.
> >>>>>>>>>> Tech-X Corp. Research Scientist
> >>>>>>>>>> 5621 Arapahoe Ave. Suite A
> >>>>>>>>>> Boulder, CO 80303
> >>>>>>>>>> 303-996-2034
> >>>>>>>>>> =========================
> >>>>>>>>>>
> >>>>>>>>>
> >>>>>>>>
> >>>>>>>>
> >>>>>>>
> >>>>>>
> >>>>>>
> >>>>>>
> >>>>>>
> >>>>>
> >>>>
> >>>>
> >>>>
> >>>>
> >>
> >>
> >
>
>
More information about the petsc-dev
mailing list