32 throw ProSHADE_exception (
"Attempted allocating Wigner D matrices before\n : allocating E matrices memory.",
"EW00024", __FILE__, __LINE__, __func__,
"The E matrices and Wigner matrices both require to know\n : the bandwidth of the comparison (which may differ from the\n : object bandwidth). This is set when allocating E matrices\n : and therefore if it is 0 now, E matrices were not yet\n : allocated." );
40 for ( proshade_unsign bandIter = 0; bandIter < this->
maxCompBand; bandIter++ )
43 this->
wignerMatrices[bandIter] =
new proshade_complex* [(bandIter * 2) + 1];
47 for ( proshade_unsign order1Iter = 0; order1Iter < ( (bandIter * 2) + 1 ); order1Iter++ )
49 this->
wignerMatrices[bandIter][order1Iter] =
new proshade_complex [(bandIter * 2) + 1];
75 void ProSHADE_internal_wigner::allocateWignerWorkspace ( proshade_double*& matIn, proshade_double*& matOut, proshade_double*& sqrts, proshade_double*& workspace, proshade_double*& alphaExponentReal, proshade_double*& alphaExponentImag, proshade_double*& gammaExponentReal, proshade_double*& gammaExponentImag, proshade_double*& trigs, proshade_unsign compBand )
78 matIn =
new proshade_double[
static_cast<proshade_unsign
> ( 4 * pow( compBand, 2.0 ) ) - 4 * compBand + 1];
80 matOut =
new proshade_double[
static_cast<proshade_unsign
> ( 4 * pow( compBand, 2.0 ) ) - 4 * compBand + 1];
81 sqrts =
new proshade_double[
static_cast<proshade_unsign
> ( 2 * compBand )];
82 workspace =
new proshade_double[
static_cast<proshade_unsign
> ( 4 * pow( compBand, 2.0 ) )];
83 alphaExponentReal =
new proshade_double[
static_cast<proshade_unsign
> ( 2 * compBand - 1 )];
84 alphaExponentImag =
new proshade_double[
static_cast<proshade_unsign
> ( 2 * compBand - 1 )];
85 gammaExponentReal =
new proshade_double[
static_cast<proshade_unsign
> ( 2 * compBand - 1 )];
86 gammaExponentImag =
new proshade_double[
static_cast<proshade_unsign
> ( 2 * compBand - 1 )];
87 trigs =
new proshade_double[2];
116 void ProSHADE_internal_wigner::releaseWignerWorkspace ( proshade_double*& matIn, proshade_double*& matOut, proshade_double*& sqrts, proshade_double*& workspace, proshade_double*& alphaExponentReal, proshade_double*& alphaExponentImag, proshade_double*& gammaExponentReal, proshade_double*& gammaExponentImag, proshade_double*& trigs )
119 if ( matIn !=
nullptr ) {
delete[] matIn; }
120 if ( matOut !=
nullptr ) {
delete[] matOut; }
121 if ( sqrts !=
nullptr ) {
delete[] sqrts; }
122 if ( workspace !=
nullptr ) {
delete[] workspace; }
123 if ( trigs !=
nullptr ) {
delete[] trigs; }
124 if ( alphaExponentReal !=
nullptr ) {
delete[] alphaExponentReal; }
125 if ( alphaExponentImag !=
nullptr ) {
delete[] alphaExponentImag; }
126 if ( gammaExponentReal !=
nullptr ) {
delete[] gammaExponentReal; }
127 if ( gammaExponentImag !=
nullptr ) {
delete[] gammaExponentImag; }
147 void ProSHADE_internal_wigner::prepareTrigsSqrtsAndExponents ( proshade_double* sqrts, proshade_double* alphaExponentReal, proshade_double* alphaExponentImag, proshade_double* gammaExponentReal, proshade_double* gammaExponentImag, proshade_double* trigs, proshade_unsign compBand, proshade_double angAlpha, proshade_double angBeta, proshade_double angGamma )
150 for ( proshade_unsign iter = 0; iter < ( 2 * compBand ); iter++ )
152 sqrts[iter] =
static_cast<proshade_double
> ( sqrt (
static_cast<proshade_double
> ( iter ) ) );
156 trigs[0] =
static_cast<proshade_double
> ( cos ( 0.5 * -angBeta ) );
157 trigs[1] =
static_cast<proshade_double
> ( sin ( 0.5 * -angBeta ) );
160 genExp (
static_cast< int > ( compBand ), angAlpha, alphaExponentReal, alphaExponentImag );
161 genExp (
static_cast< int > ( compBand ), angGamma, gammaExponentReal, gammaExponentImag );
187 void ProSHADE_internal_wigner::computeWignerMatrices (
ProSHADE_settings* settings,
ProSHADE_internal_data::ProSHADE_data* obj, proshade_double* alphaExponentReal, proshade_double* alphaExponentImag, proshade_double* gammaExponentReal, proshade_double* gammaExponentImag, proshade_double* matIn, proshade_double* matOut, proshade_double* trigs, proshade_double* sqrts, proshade_double* workspace )
193 proshade_double *expARStart, *expAIStart, *expGRStart, *expGIStart;
194 proshade_double Dij, eARi, eAIi, eGRj, eGIj, iSign, rSign;
195 proshade_complex hlpVal;
196 proshade_unsign noOrders, arrConvIter;
197 for ( proshade_unsign bandIter = 0; bandIter < obj->
getComparisonBand(); bandIter++ )
200 noOrders = 2 * bandIter + 1;
210 wignerdmat (
static_cast< int > ( bandIter ), matIn, matOut, trigs, sqrts, workspace );
213 for ( proshade_unsign d1Iter = 0; d1Iter < noOrders; d1Iter++ )
215 eARi = expARStart[d1Iter];
216 eAIi = expAIStart[d1Iter];
218 for ( proshade_unsign d2Iter = 0; d2Iter < noOrders; d2Iter++ )
220 Dij = matOut[arrConvIter];
221 eGRj = expGRStart[d2Iter];
222 eGIj = expGIStart[d2Iter];
224 hlpVal[0] = ( Dij * eGRj * eARi - Dij * eGIj * eAIi ) * rSign;
225 hlpVal[1] = ( Dij * eGRj * eAIi + Dij * eGIj * eARi ) * iSign;
235 memcpy ( matIn, matOut,
sizeof ( proshade_double ) * ( noOrders * noOrders ) );
261 proshade_double *matIn, *matOut, *sqrts, *workspace, *alphaExponentReal, *alphaExponentImag, *gammaExponentReal, *gammaExponentImag, *trigs;
275 computeWignerMatrices ( settings, obj, alphaExponentReal, alphaExponentImag, gammaExponentReal, gammaExponentImag,
276 matIn, matOut, trigs, sqrts, workspace );
280 gammaExponentReal, gammaExponentImag, trigs );