1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
|
// ************************************************************************************************
//
// BornAgain: simulate and fit reflection and scattering
//
//! @file Resample/Interparticle/SSCAStrategy.cpp
//! @brief Implements class SSCAStrategy.
//!
//! @homepage http://www.bornagainproject.org
//! @license GNU General Public License v3 or higher (see COPYING)
//! @copyright Forschungszentrum Jülich GmbH 2018
//! @authors Scientific Computing Group at MLZ (see CITATION, AUTHORS)
//
// ************************************************************************************************
#include "Resample/Interparticle/SSCAStrategy.h"
#include "Resample/Coherence/CoheringSubparticles.h"
#include "Resample/Element/DiffuseElement.h"
#include "Sample/Aggregate/InterferenceRadialParacrystal.h"
namespace {
double meanRadius(const OwningVector<const CoheringSubparticles>& weighted_formfactors)
{
double result = 0.0;
for (const auto& ffw : weighted_formfactors)
result += ffw->relativeAbundance() * ffw->radialExtension();
return result;
}
} // namespace
SSCAStrategy::SSCAStrategy(const OwningVector<const CoheringSubparticles>& weighted_formfactors,
const InterferenceRadialParacrystal* iff, SimulationOptions options,
bool polarized, double kappa)
: IInterparticleStrategy(weighted_formfactors, options, polarized)
, m_iff(iff->clone())
, m_kappa(kappa)
, m_mean_radius(meanRadius(weighted_formfactors))
{
}
//! Returns the total scattering intensity for given kf and
//! for one particle layout (implied by the given particle form factors).
//! This is the scalar version
double SSCAStrategy::scalarCalculation(const DiffuseElement& ele) const
{
const double qp = ele.meanQ().magxy();
double diffuse_intensity = 0.0;
complex_t ff_orig = 0., ff_conj = 0.; // original and conjugated mean formfactor
for (const auto& ffw : m_weighted_formfactors) {
complex_t ff = ffw->summedFF(ele);
double fraction = ffw->relativeAbundance();
diffuse_intensity += fraction * std::norm(ff);
double radial_extension = ffw->radialExtension();
complex_t prefac =
ffw->relativeAbundance() * calculatePositionOffsetPhase(qp, radial_extension);
ff_orig += prefac * ff;
ff_conj += prefac * std::conj(ff);
}
const complex_t mean_ff_norm = ff_orig * ff_conj;
const complex_t p2kappa = getCharacteristicSizeCoupling(qp, m_weighted_formfactors);
const complex_t omega = m_iff->FTPDF(qp);
const double iff = 2.0 * (mean_ff_norm * omega / (1.0 - p2kappa * omega)).real();
const double dw_factor = m_iff->DWfactor(ele.meanQ());
return diffuse_intensity + dw_factor * iff;
}
//! This is the polarized version
double SSCAStrategy::polarizedCalculation(const DiffuseElement& ele) const
{
const double qp = ele.meanQ().magxy();
SpinMatrix diffuse_matrix;
const SpinMatrix& polarizer = ele.polarizer();
const SpinMatrix& analyzer = ele.analyzer();
SpinMatrix ff_orig;
SpinMatrix ff_conj;
for (const auto& ffw : m_weighted_formfactors) {
const SpinMatrix ff = ffw->summedPolFF(ele);
const double fraction = ffw->relativeAbundance();
diffuse_matrix += fraction * (ff * polarizer * ff.adjoint());
const double radial_extension = ffw->radialExtension();
const complex_t prefac =
ffw->relativeAbundance() * calculatePositionOffsetPhase(qp, radial_extension);
ff_orig += prefac * ff;
ff_conj += prefac * ff.adjoint();
}
const complex_t p2kappa = getCharacteristicSizeCoupling(qp, m_weighted_formfactors);
const complex_t omega = m_iff->FTPDF(qp);
const SpinMatrix interference_matrix =
(2.0 * omega / (1.0 - p2kappa * omega)) * analyzer * ff_orig * polarizer * ff_conj;
const SpinMatrix diffuse_matrix2 = analyzer * diffuse_matrix;
const double interference_trace = std::abs(interference_matrix.trace());
const double diffuse_trace = std::abs(diffuse_matrix2.trace());
const double dw_factor = m_iff->DWfactor(ele.meanQ());
return diffuse_trace + dw_factor * interference_trace;
}
complex_t SSCAStrategy::getCharacteristicSizeCoupling(
double qp, const OwningVector<const CoheringSubparticles>& ff_wrappers) const
{
complex_t result = 0;
for (const auto& ffw : ff_wrappers)
result += ffw->relativeAbundance()
* calculatePositionOffsetPhase(2.0 * qp, ffw->radialExtension());
return result;
}
complex_t SSCAStrategy::calculatePositionOffsetPhase(double qp, double radial_extension) const
{
return exp_I(m_kappa * qp * (radial_extension - m_mean_radius));
}
|
