Menu
▾
▴
cppqed-support — Bugs, questions, feature requests, etc.
You can subscribe to this list here.
| 2012 |
Jan
|
Feb
(5) |
Mar
(5) |
Apr
(10) |
May
(5) |
Jun
|
Jul
|
Aug
(17) |
Sep
(6) |
Oct
(1) |
Nov
|
Dec
|
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 2013 |
Jan
(1) |
Feb
(4) |
Mar
(2) |
Apr
|
May
|
Jun
(1) |
Jul
|
Aug
(3) |
Sep
|
Oct
(1) |
Nov
(18) |
Dec
(10) |
| 2014 |
Jan
(12) |
Feb
|
Mar
(2) |
Apr
(11) |
May
(1) |
Jun
|
Jul
|
Aug
|
Sep
|
Oct
(3) |
Nov
|
Dec
|
|
From: Raimar S. <rai...@ui...> - 2012-04-23 11:49:50
|
Dear András, thank you for adjusting and merging the fix. On Saturday 21 April 2012 10:28:58 you wrote: > Great! Thanks a lot for tracking this down. But this means that the first > timestep after a continuation was always erroneous so far (that is, if the > system time development had *Exact* part)? And this eventually produced an > overflow when the trajectory was continued after long enough a time? The stack when psi first gets corrupted is something like (left out some indirections): quantumtrajectory::evolve -> MCWF_Trajectory<RANK>::step -> MCWF_Trajectory<RANK>::coherentTimeDevelopment -> Evolved<A>::step -> evolved::details::apply -> gsl_odeiv_evolve_apply -> MCWF_Trajectory<RANK>::derivs -> Hamiltonian::addContribution -> BinarySystem::addContribution -> Hamiltonian<2,ONE_TIME>::addContribution Here, after the call to addContribution of the binary Hamiltonian (the frees are ok in the first step, why?), the maximum absolute value in dpsidt is 3.5e+130. I would think this is because the transformed non-hermitian Hamiltonian has some huge elements as we evaluate U(t) for large t. These large values carry on to psi as the odeiv stepper takes the huge derivatives and updates the state (max is 4e+126). The next time derivs is called, the maximum absolute value of dpsidt is already on the order of 3e+257 and psi gets updated to something on the order of 3e+257. Then first dpsidt and afterwards psi overflows. What I don't understand is why derivs gets called three times with different values of t (t=30, t=30.00014 and t=30.00021) before dpsidt overflows, but each time tIntPic0_=0. Shouldn't tIntPic0_ get updated? And the question remains if continuation always was erroneous. I would argue that as long as the stepper makes it to the first update of tIntPic0_ without overflow, everything should be fine except if rounding errors skew the results. Physically large values of t should be legal, the problems arise if we cannot represent the transformation numerically accurate. Best regards Raimar |
|
From: Andras V. <and...@ui...> - 2012-04-21 08:29:25
|
Dear Raimar, Great! Thanks a lot for tracking this down. But this means that the first timestep after a continuation was always erroneous so far (that is, if the system time development had *Exact* part)? And this eventually produced an overflow when the trajectory was continued after long enough a time? A little correction is yet necessary, however: if you look at * MCWF_Trajectory<RANK>::coherentTimeDevelopment*, you can see that *tIntPic0_ * is expected to equal the beginning of the present timestep only if *ex_*is true. Otherwise, it is expected to remain zero. (You can find a tentative explanation for this here<http://cppqed.sourceforge.net/structure.html#hamiltonian> .) This condition has to be included here as well. I have merged, corrected and pushed revision #212. Could you please have a look at * MCWF_Trajectory.tcc:128–130* once more and see whether you agree? (You can also remove the raimar/continue_fix branch). Thank again, and best regards, András Dr. Andras Vukics Institute for Theoretical Physics University of Innsbruck |
|
From: Raimar S. <rai...@ui...> - 2012-04-20 13:54:31
|
Dear András, I think I found a fix to this problem, pleas have a look at the raimar/continue_fix branch. To reproduce the problem and test the fix it is sufficient to do 1particle1mode --o try2.d --T 30 --dc 0 --Dt 1 1particle1mode --o try2.d --T 32 --dc 0 --Dt 1 Best regards, Raimar On Thursday 19 April 2012 10:48:31 Andras Vukics wrote: > Dear Raimar, > > Actually, no, because the timestep displayed in the 2nd column is dtDid and > not dtTry, while the timestep read in from the .sv file is dtTry. This was > my confusion, too. > I have pushed a new revision recently, with a bit clearer output in the > case of trajectory continuation to remove this confusion. > > So, I still (or, again) have no idea why you get these infinites. > > Keep in touch, > András > > > > > On Thu, Apr 19, 2012 at 10:30 AM, Raimar Sandner > > <rai...@ui...>wrote: > > Dear András, > > > > I want to look into this problem a little bit more. You closed the bug > > about > > the erroneous parsing of the timestep as invalid, but it seems to be a > > valid > > bug. The timestep is read in as 1 which is clearly wrong or am I mistaken? > > > > Best regards > > Raimar > > > > On Friday 30 March 2012 00:03:46 Andras Vukics wrote: > > > I could already detect the error at T=30 (after several continuations). > > > There is also a hint of some bug, cf. > > > > http://sourceforge.net/tracker/?func=detail&aid=3512988&group_id=187775&at > > id> > > > =922653 > > > On Thu, Mar 29, 2012 at 9:16 PM, Andras Vukics > > > > <and...@ui...>wrote: > > > > Dear Raimar, > > > > I've started to experiment, and try to reproduce your error. > > > > Just a quick question: did you try in debug mode? > > > > András > > > > -------------------------------------------------------------------------- > > ---- For Developers, A Lot Can Happen In A Second. > > Boundary is the first to Know...and Tell You. > > Monitor Your Applications in Ultra-Fine Resolution. Try it FREE! > > http://p.sf.net/sfu/Boundary-d2dvs2 > > _______________________________________________ > > Cppqed-support mailing list > > Cpp...@li... > > https://lists.sourceforge.net/lists/listinfo/cppqed-support |
|
From: Andras V. <and...@ui...> - 2012-04-19 08:49:02
|
Dear Raimar, Actually, no, because the timestep displayed in the 2nd column is dtDid and not dtTry, while the timestep read in from the .sv file is dtTry. This was my confusion, too. I have pushed a new revision recently, with a bit clearer output in the case of trajectory continuation to remove this confusion. So, I still (or, again) have no idea why you get these infinites. Keep in touch, András On Thu, Apr 19, 2012 at 10:30 AM, Raimar Sandner <rai...@ui...>wrote: > Dear András, > > I want to look into this problem a little bit more. You closed the bug > about > the erroneous parsing of the timestep as invalid, but it seems to be a > valid > bug. The timestep is read in as 1 which is clearly wrong or am I mistaken? > > Best regards > Raimar > > On Friday 30 March 2012 00:03:46 Andras Vukics wrote: > > I could already detect the error at T=30 (after several continuations). > > There is also a hint of some bug, cf. > > > http://sourceforge.net/tracker/?func=detail&aid=3512988&group_id=187775&atid > > =922653 > > On Thu, Mar 29, 2012 at 9:16 PM, Andras Vukics > <and...@ui...>wrote: > > > Dear Raimar, > > > I've started to experiment, and try to reproduce your error. > > > Just a quick question: did you try in debug mode? > > > András > > > > > > > > > > > > > ------------------------------------------------------------------------------ > For Developers, A Lot Can Happen In A Second. > Boundary is the first to Know...and Tell You. > Monitor Your Applications in Ultra-Fine Resolution. Try it FREE! > http://p.sf.net/sfu/Boundary-d2dvs2 > _______________________________________________ > Cppqed-support mailing list > Cpp...@li... > https://lists.sourceforge.net/lists/listinfo/cppqed-support > |
|
From: Raimar S. <rai...@ui...> - 2012-04-19 08:31:48
|
Dear András, I want to look into this problem a little bit more. You closed the bug about the erroneous parsing of the timestep as invalid, but it seems to be a valid bug. The timestep is read in as 1 which is clearly wrong or am I mistaken? Best regards Raimar On Friday 30 March 2012 00:03:46 Andras Vukics wrote: > I could already detect the error at T=30 (after several continuations). > There is also a hint of some bug, cf. > http://sourceforge.net/tracker/?func=detail&aid=3512988&group_id=187775&atid > =922653 > On Thu, Mar 29, 2012 at 9:16 PM, Andras Vukics <and...@ui...>wrote: > > Dear Raimar, > > I've started to experiment, and try to reproduce your error. > > Just a quick question: did you try in debug mode? > > András > > > > > > |
|
From: Raimar S. <rai...@ui...> - 2012-04-18 13:28:47
|
Dear András,
I have made some progress trying to compile the framework with Intel's compiler.
The reason why I try is that I have an account at the university's leo3 HPC with Intel
processors, and I hope for increased performance.
In debug mode, the testsuitePhysics is passed. In release mode, testsuitePhysics is passed
except tests 4 and 64, where I get SIGBUS or SIGSEGV depending on the intel compiler version.
This happens in cpputils::details::concatenateHelper upon return std::copy(t.begin(),t.end(),iter);
Unfortunately my own script also crashes with SIGBUS.
Preliminary tests show that the scripts which succeed are indeed faster, though I have to investigate
this more systematically.
I will also try to find the exact compiler flags which cause the crashes.
Best regards
Raimar
P.S. let me know if you want to experiment with this on teazer, I have a version of gsl compiled
with the intel compiler and I can send you the user-config.jam which points to this library.
The changes that were necessary to compile the code are in the branch cppqed/raimar/intel.
Here is a summary:
1. namespace problems which icpc has for some reasons
========================================================================
2. Undefined references to some static const members, fixed by e.g.
=== modified file 'utils/include/impl/BlitzArraySliceIterator.tcc'
--- utils/include/impl/BlitzArraySliceIterator.tcc 2011-09-19 11:29:15 +0000
+++ utils/include/impl/BlitzArraySliceIterator.tcc 2012-04-18 11:49:57 +0000
@@ -205,6 +205,9 @@
template<int RANK, typename V>
typename IndexerBase<RANK,V>::Idx IndexerBase<RANK,V>::cache_;
+template<int RANK, typename V>
+const typename IndexerBase<RANK,V>::Helper IndexerBase<RANK,V>::helper_ = IndexerBase<RANK,V>::Helper();
+
========================================================================
3. Here I am not sure why gcc compiles the code, as the constructor which is called is missing:
=== modified file 'elements/frees/MultiLevel.h'
--- elements/frees/MultiLevel.h 2011-08-11 13:21:25 +0000
+++ elements/frees/MultiLevel.h 2012-04-18 11:49:57 +0000
@@ -316,7 +316,7 @@
// here, the compiler cannot deduce NL
// anymore
const VP& etas, const VL& gammas,
- const Averaged& averaged=multilevel::DiagonalDO(NL))
+ const Averaged& averaged=multilevel::DiagonalDO("",NL))
{
return RETURN_type(new PumpedLossyMultiLevelSch<NL,VP,VL,Averaged>(deltas,etas,gammas,averaged));
}
@@ -325,7 +325,7 @@
template<int NL, typename VP, typename VL, typename Averaged>
inline
RETURN_type
-makePumpedLossyMultiLevelSch(const multilevel::ParsPumpedLossy<NL,VP,VL>& p, const Averaged& averaged=multilevel::DiagonalDO(NL))
+makePumpedLossyMultiLevelSch(const multilevel::ParsPumpedLossy<NL,VP,VL>& p, const Averaged& averaged=multilevel::DiagonalDO("",NL))
{
return RETURN_type(new PumpedLossyMultiLevelSch<NL,VP,VL,Averaged>(p.deltas,p.etas,p.gammas,averaged));
}
=========================================================================
4. Here I have no clue what is going on, adding the public keyword to a struct should not make any difference.
But icpc complains with "... is inaccessible" without it.
=== modified file 'elements/interactions/impl/MLJC.tcc'
--- elements/interactions/impl/MLJC.tcc 2011-10-07 23:10:35 +0000
+++ elements/interactions/impl/MLJC.tcc 2012-04-18 11:49:57 +0000
@@ -23,6 +23,7 @@
template<int NL, typename VC> template<int N1, int N2>
struct Base<NL,VC>::ModeDynamics : public tmptools::pair_c<N1,N2>
{
+public:
ModeDynamics(const MultiLevelBase<NL>& ml, const ModeBase& mode, const dcomp& g)
: a (-g*mode::aop(&mode)),
adagger((g*mode::aop(&mode)).dagger())
|
|
From: Andras V. <and...@ui...> - 2012-03-29 22:04:07
|
I could already detect the error at T=30 (after several continuations). There is also a hint of some bug, cf. http://sourceforge.net/tracker/?func=detail&aid=3512988&group_id=187775&atid=922653 On Thu, Mar 29, 2012 at 9:16 PM, Andras Vukics <and...@ui...>wrote: > Dear Raimar, > I've started to experiment, and try to reproduce your error. > Just a quick question: did you try in debug mode? > András > > > > On Wed, Mar 21, 2012 at 6:59 PM, Raimar Sandner <rai...@ui... > > wrote: > >> Dear András, >> >> I tried to continue a MCWF simulation by increasing --T and leaving the >> output >> file and .sv file in place. The simulation aborted with an >> InfiniteDetectedException. >> >> Indeed: I looked into the average function of the mode, the first >> displayed >> step is fine, but in the second display step this is the >> LazyDensityOperator >> which is passed to the average method (std::cout of the corresponding >> DensityOperatorLow): >> >> (0,4) x (0,4) >> [ (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) >> (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) >> (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) >> (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) >> (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) ] >> >> It is the same for the other frees of the system. Unfortunately I >> couldn't find >> the reason. >> >> This only happens for longer trajectories, i.e. everything is fine if I >> continue a short trajectory with T=15, but I get the exception with T=120. >> >> For the case that you don't know what is going on right away, I have >> attached >> the output file and .sv file as a test case. You would need to checkout >> the >> raimar/private bazaar branch for the script and the interactions I use >> therein. The command line is >> >> 2ParticlesRingCavityFull --etaSin 9 --minitSin 0.1 --modeCavSin Sin \ >> --kappaCos 5 --pinit1 "(0.5 0 0.09 0)" --pinit2 "(-0.5 0 0.09 0)" \ >> --UnotSin -5 --Dt 1 --cutoffCos 3 --deltaCSin -20 --deltaCCos -20 \ >> --evol single --fin1 5 --fin2 5 --precision 6 --dc 0 --minitCos 0.1 \ >> --cutoffSin 5 --modeCavCos Cos --UnotCos -5 --kappaSin 5 --seed 1001 \ >> --o 2ParticlesRingCavityFull.out.1001 --T 121 >> >> Best regards, >> Raimar >> >> ------------------------------------------------------------------------------ >> This SF email is sponsosred by: >> Try Windows Azure free for 90 days Click Here >> http://p.sf.net/sfu/sfd2d-msazure >> _______________________________________________ >> Cppqed-support mailing list >> Cpp...@li... >> https://lists.sourceforge.net/lists/listinfo/cppqed-support >> >> > |
|
From: Andras V. <and...@ui...> - 2012-03-29 19:17:23
|
Dear Raimar, I've started to experiment, and try to reproduce your error. Just a quick question: did you try in debug mode? András On Wed, Mar 21, 2012 at 6:59 PM, Raimar Sandner <rai...@ui...>wrote: > Dear András, > > I tried to continue a MCWF simulation by increasing --T and leaving the > output > file and .sv file in place. The simulation aborted with an > InfiniteDetectedException. > > Indeed: I looked into the average function of the mode, the first displayed > step is fine, but in the second display step this is the > LazyDensityOperator > which is passed to the average method (std::cout of the corresponding > DensityOperatorLow): > > (0,4) x (0,4) > [ (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) > (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) > (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) > (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) > (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) ] > > It is the same for the other frees of the system. Unfortunately I couldn't > find > the reason. > > This only happens for longer trajectories, i.e. everything is fine if I > continue a short trajectory with T=15, but I get the exception with T=120. > > For the case that you don't know what is going on right away, I have > attached > the output file and .sv file as a test case. You would need to checkout the > raimar/private bazaar branch for the script and the interactions I use > therein. The command line is > > 2ParticlesRingCavityFull --etaSin 9 --minitSin 0.1 --modeCavSin Sin \ > --kappaCos 5 --pinit1 "(0.5 0 0.09 0)" --pinit2 "(-0.5 0 0.09 0)" \ > --UnotSin -5 --Dt 1 --cutoffCos 3 --deltaCSin -20 --deltaCCos -20 \ > --evol single --fin1 5 --fin2 5 --precision 6 --dc 0 --minitCos 0.1 \ > --cutoffSin 5 --modeCavCos Cos --UnotCos -5 --kappaSin 5 --seed 1001 \ > --o 2ParticlesRingCavityFull.out.1001 --T 121 > > Best regards, > Raimar > > ------------------------------------------------------------------------------ > This SF email is sponsosred by: > Try Windows Azure free for 90 days Click Here > http://p.sf.net/sfu/sfd2d-msazure > _______________________________________________ > Cppqed-support mailing list > Cpp...@li... > https://lists.sourceforge.net/lists/listinfo/cppqed-support > > |
|
From: Andras V. <and...@ui...> - 2012-03-23 09:14:34
|
Dear Raimar, It is the set of retained index *positions* which is known at compile time. The set of dummy index *positions* is indeed just the complement of this. But the set of dummy indices (i.e., their *values*) becomes available only at runtime. This is how I meant it. Then, the *values* of the retained indices is of course again a runtime-information, but this is used only when we actually address certain elements of the *slice*. The slice itself, however, is determined by the set of retained index positions and the set of dummy index values, that is, a piece of compile-time and a piece of runtime information. I hope it is clearer now. Best regards, András On Fri, Mar 23, 2012 at 9:55 AM, Raimar Sandner <rai...@ui...>wrote: > Dear András, > > just a quick question: > > In section 4 of your recent paper you state that the set of retained index > positions is an information available at compile time, while the set of > dummy > indices is an information becoming available only at runtime. > > Why is that? This is not clear to me, as the layout of the system is known > at > compile time, and I would think the set of dummy indices is just the > complement of the set of the retained indices in the set of all indices. > > Best regards > Raimar > > > ------------------------------------------------------------------------------ > This SF email is sponsosred by: > Try Windows Azure free for 90 days Click Here > http://p.sf.net/sfu/sfd2d-msazure > _______________________________________________ > Cppqed-support mailing list > Cpp...@li... > https://lists.sourceforge.net/lists/listinfo/cppqed-support > |
|
From: Raimar S. <rai...@ui...> - 2012-03-23 08:56:59
|
Dear András, just a quick question: In section 4 of your recent paper you state that the set of retained index positions is an information available at compile time, while the set of dummy indices is an information becoming available only at runtime. Why is that? This is not clear to me, as the layout of the system is known at compile time, and I would think the set of dummy indices is just the complement of the set of the retained indices in the set of all indices. Best regards Raimar |
|
From: Raimar S. <rai...@ui...> - 2012-03-21 18:00:13
|
Dear András, I tried to continue a MCWF simulation by increasing --T and leaving the output file and .sv file in place. The simulation aborted with an InfiniteDetectedException. Indeed: I looked into the average function of the mode, the first displayed step is fine, but in the second display step this is the LazyDensityOperator which is passed to the average method (std::cout of the corresponding DensityOperatorLow): (0,4) x (0,4) [ (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) (-nan,-nan) ] It is the same for the other frees of the system. Unfortunately I couldn't find the reason. This only happens for longer trajectories, i.e. everything is fine if I continue a short trajectory with T=15, but I get the exception with T=120. For the case that you don't know what is going on right away, I have attached the output file and .sv file as a test case. You would need to checkout the raimar/private bazaar branch for the script and the interactions I use therein. The command line is 2ParticlesRingCavityFull --etaSin 9 --minitSin 0.1 --modeCavSin Sin \ --kappaCos 5 --pinit1 "(0.5 0 0.09 0)" --pinit2 "(-0.5 0 0.09 0)" \ --UnotSin -5 --Dt 1 --cutoffCos 3 --deltaCSin -20 --deltaCCos -20 \ --evol single --fin1 5 --fin2 5 --precision 6 --dc 0 --minitCos 0.1 \ --cutoffSin 5 --modeCavCos Cos --UnotCos -5 --kappaSin 5 --seed 1001 \ --o 2ParticlesRingCavityFull.out.1001 --T 121 Best regards, Raimar |
|
From: Andras V. <and...@ui...> - 2012-02-27 20:39:26
|
Dear "members of the C++QED community", I'm pleased to announce that my paper on C++QEDv2 has been published: http://dx.doi.org/10.1016/j.cpc.2012.02.004 Best regards, András |
|
From: Andras V. <and...@ui...> - 2012-02-15 14:28:52
|
Hi Raimar, thanks for your reply. So you are saying blitzplusplus::vfmsi::fullRange > could > always only give me the two-particle (left or right) state vector from my > DensityOperatorLow, this is why I have to use > blitzplusplus::basi::fullRange. > That's right. On a related note, would it be possible that, given a > blitzplusplus::basi::Iterator with two or more retained index positions, to > slice it even further? This would somehow be related to the concept that a > slice of a multi-array is again a multi-array and can be sliced. > Yes, a sliced multi-array (e.g. a dereferenced slice iterator) can be further sliced (can be used to initialize a slice iterator). This actually does occur in ModeCorrelations::average, which in general gets a sliced LazyDensityOperator, and slices it further when calculating the averages for the individual Modes relying on mode::AveragedQuadratures (cf. ModeCorrelations.cc:38-39). You should take ModeCorrelations as an example, because also with the particle correlations, some functionalities could (perhaps) be delegated to particle::Averaged. So in principle the slicing would occur in > the fft function, not in the average function anymore. How would such a fft > function with specialized implementation determine the type of the passed > LazyDensityOperator, by dynamic_cast? > Yes, I think this is a paradigmatic example of dynamic_cast. The function should have roughly this signature: template<int RANK, int IDX> const LazyDensityOperator<RANK> ffTransformLazyDensityOperator(const LazyDensityOperator<RANK>&, fft::Direction, mpl::int_<IDX>); IDX is the index position on which the transformation should occur. And in the density-operator case, the transformation occurs on the index position RANK+IDX as well, only with opposit Copying is of course still necessary, because we do not want to touch the original one (in fact, we cannot). But in case of a state vector, only this needs to be copied. Wouldn't boost::mpl::range_c<int,0,1>, boost::mpl::range_c<int,1,2> etc > work, > although of course a bit awkward? > They should indeed work, but it's an overkill :) . Best regards, András |
|
From: Raimar S. <rai...@ui...> - 2012-02-15 11:48:45
|
Dear András, thanks for your reply. So you are saying blitzplusplus::vfmsi::fullRange could always only give me the two-particle (left or right) state vector from my DensityOperatorLow, this is why I have to use blitzplusplus::basi::fullRange. Actually, your code snippet is pretty much what I also had in mind. On a related note, would it be possible that, given a blitzplusplus::basi::Iterator with two or more retained index positions, to slice it even further? This would somehow be related to the concept that a slice of a multi-array is again a multi-array and can be sliced. On Wednesday 15 February 2012 11:53:02 Andras Vukics wrote: > This is actually rather wasteful, because even in the case when the > underlying data is just a state vector, it is treated in the same way as a > full density operator. You could maybe think about a nice general solution > to this problem (I have just added a tracker #3487814 on this issue). Yes I have noticed the tracker :). So in principle the slicing would occur in the fft function, not in the average function anymore. How would such a fft function with specialized implementation determine the type of the passed LazyDensityOperator, by dynamic_cast? > // left: > for_each(fullRange(rhotemp,Vector<0>()),bind(&ffTransform,_1,FFTDIR_KX)); > for_each(fullRange(rhotemp,Vector<1>()),bind(&ffTransform,_1,FFTDIR_KX)); > // right: > for_each(fullRange(rhotemp,Vector<2>()),bind(&ffTransform,_1,FFTDIR_XK)); > for_each(fullRange(rhotemp,Vector<3>()),bind(&ffTransform,_1,FFTDIR_XK)); > > Where Vector stands for tmptools::Vector. (This is not the only class that > works with the slice iterators, because e.g. boost::mpl::range_c also > works, but here you need this.) Wouldn't boost::mpl::range_c<int,0,1>, boost::mpl::range_c<int,1,2> etc work, although of course a bit awkward? Best regards Raimar |
|
From: Andras V. <and...@ui...> - 2012-02-15 10:53:33
|
Hi Raimar, In this case, your *LazyDensityOperator* will be binary only, which means that it can be indexed with 4 indices (2 instances of *TTD_IDXTINY(2)* ). You have to copy the whole thing into a * quantumdata::Types<2>::DensityOperatorLow* . This is actually rather wasteful, because even in the case when the underlying data is just a state vector, it is treated in the same way as a full density operator. You could maybe think about a nice general solution to this problem (I have just added a tracker #3487814 on this issue). Anyway, once you have your *DensityOperatorLow*, I think you have to apply FFT on *all of its indices individually,* the FFT's sign depending on whether it is a ‘left’ or ‘right’ index. So, * blitzplusplus::vfmsi::fullRange* will not be of much use in this case, and you have to apply the more general *blitzplusplus::basi::fullRange*iteration. I think the code is something pretty much like this: // left: for_each(fullRange(rhotemp,Vector<0>()),bind(&ffTransform,_1,FFTDIR_KX)); for_each(fullRange(rhotemp,Vector<1>()),bind(&ffTransform,_1,FFTDIR_KX)); // right: for_each(fullRange(rhotemp,Vector<2>()),bind(&ffTransform,_1,FFTDIR_XK)); for_each(fullRange(rhotemp,Vector<3>()),bind(&ffTransform,_1,FFTDIR_XK)); Where *Vector* stands for *tmptools::Vector*. (This is not the only class that works with the slice iterators, because e.g. *boost::mpl::range_c*also works, but here you need this.) Keep in touch. Best regards, András |
|
From: Raimar S. <rai...@ui...> - 2012-02-10 17:26:30
|
Dear András, I want to implement an Averaged class which calculates the position correlation of two particles. For this I have to calculate the Fourier transform of a *LazyDensityOperator* with rank 4. In *Particle.cc* you copy the *LazyDensityOperator* into a temporary *DensityOperatorLow*. I assume this is necessary because the gsl FFT doesn't work directly with a *LazyDensityOperator*? >From my understanding for_each(fullRange(rhotemp, Left()),bind(&ffTransform,_1,FFTDIR_KX)); for_each(fullRange(rhotemp,Right()),bind(&ffTransform,_1,FFTDIR_XK)); calculates the FFT with +- ikx on the rows/columns of the density matrix, respectively. I think in my case I would have to created slices which correspond to first particle left, second particle left, first particle right and second particle right and then FFT all these slices. Can you give me a hint how this could be implemented? Should I slice with *Vector<0>()* ... *Vector<3>* directly as I know my density operator will have rank 4, or is there a more generic way of doing this similar to Left() Right()? Thanks and very best regards Raimar |
2 messages has been excluded from this view by a project administrator.
