From 5166ec0b90953300744204c690187337fdc8d09f Mon Sep 17 00:00:00 2001 From: Sebastian Date: Tue, 15 Dec 2020 10:57:03 +0100 Subject: [PATCH] Fixed several issues --- .../ProgressPrint/MinimalProgressPrint.c | 2 +- .../Libraries/ProgressPrint/ProgressPrint.c | 15 ++++---- .../Logic/Routines/StatisticsRoutines.c | 35 +++++++++++++++++-- .../Logic/Routines/StatisticsRoutines.h | 29 ++++++++++----- .../ModelViewer/DxModelSceneViewModel.cs | 4 +-- .../CellReferencePositionData.cs | 2 +- 6 files changed, 66 insertions(+), 21 deletions(-) diff --git a/src/McSolver/Libraries/ProgressPrint/MinimalProgressPrint.c b/src/McSolver/Libraries/ProgressPrint/MinimalProgressPrint.c index f05329b3..1c1ce4b9 100644 --- a/src/McSolver/Libraries/ProgressPrint/MinimalProgressPrint.c +++ b/src/McSolver/Libraries/ProgressPrint/MinimalProgressPrint.c @@ -57,7 +57,7 @@ void PrintMocassinSimulationFinishInfo(SCONTEXT_PARAMETER, file_t *fstream) let moveR1 = mobilityData.EnsembleMoveR1; let moveR2 = mobilityData.EnsembleMoveR2; - let dCoef = mobilityData.DiffusionCoefficient; + let dCoef = mobilityData.DiffusionCoefficientVector; let counters = statisticsData.CounterCollection; fprintf(fstream, "P:%i:Counters(Success:"FORMAT_I64(012)"|Reject:"FORMAT_I64(012)"|Block:"FORMAT_I64(012)"|EndUnstable:"FORMAT_I64(012)"|StartUnstable:"FORMAT_I64(012)")\n", i, diff --git a/src/McSolver/Libraries/ProgressPrint/ProgressPrint.c b/src/McSolver/Libraries/ProgressPrint/ProgressPrint.c index 2d73540d..4767fffe 100644 --- a/src/McSolver/Libraries/ProgressPrint/ProgressPrint.c +++ b/src/McSolver/Libraries/ProgressPrint/ProgressPrint.c @@ -101,14 +101,17 @@ static void PrintParticleStatistics(const ParticleStatistics_t* restrict statist // Prints the passed particle mobility data to the stream static void PrintParticleMobility(const ParticleMobilityData_t* restrict data, file_t*restrict fstream) { + fprintf(fstream, MC_DEFAULT_FORMAT(MC_OUTVEC_FORMAT), "Conductivity => Components", "S m^-1", + data->ConductivityVector.A, data->ConductivityVector.B, data->ConductivityVector.C); + + fprintf(fstream, MC_DEFAULT_FORMAT(MC_OUTVEC_FORMAT), "Conductivity => Normalized (z=+1)", "S m^-1", + data->NormalizedConductivityVector.A, data->NormalizedConductivityVector.B, data->NormalizedConductivityVector.C); + fprintf(fstream, MC_DEFAULT_FORMAT(MC_OUTF64_FORMAT), "Conductivity => Field direction", "S m^-1", data->TotalConductivity); - fprintf(fstream, MC_DEFAULT_FORMAT(MC_OUTF64_FORMAT), "Conductivity => Normalized (z=+1)", "S m^-1", - data->TotalConductivityPerCharge); - - fprintf(fstream, MC_DEFAULT_FORMAT(MC_OUTVEC_FORMAT), "Conductivity => FDT components", "S m^-1", - data->NernstEinsteinConductivity.A, data->NernstEinsteinConductivity.B, data->NernstEinsteinConductivity.C); + fprintf(fstream, MC_DEFAULT_FORMAT(MC_OUTVEC_FORMAT), "Mobility => D_Sigma", "m^2 s^-1", + data->DSigmaVector.A, data->DSigmaVector.B, data->DSigmaVector.C); fprintf(fstream, MC_DEFAULT_FORMAT(MC_OUTF64_FORMAT), "Mobility => Average migration rate", "Hz", data->MigrationRate); @@ -120,7 +123,7 @@ static void PrintParticleMobility(const ParticleMobilityData_t* restrict data, f data->MobilityVector.A, data->MobilityVector.B, data->MobilityVector.C); fprintf(fstream, MC_DEFAULT_FORMAT(MC_OUTVEC_FORMAT), "Movement => Diffusion coefficient", "m^2 s^-1", - data->DiffusionCoefficient.A, data->DiffusionCoefficient.B, data->DiffusionCoefficient.C); + data->DiffusionCoefficientVector.A, data->DiffusionCoefficientVector.B, data->DiffusionCoefficientVector.C); fprintf(fstream, MC_DEFAULT_FORMAT(MC_OUTVEC_FORMAT), "Movement => Ensemble R", "m^2", data->EnsembleMoveR1.A, data->EnsembleMoveR1.B, data->EnsembleMoveR1.C); diff --git a/src/McSolver/Libraries/Simulator/Logic/Routines/StatisticsRoutines.c b/src/McSolver/Libraries/Simulator/Logic/Routines/StatisticsRoutines.c index e47671b1..28d5af41 100644 --- a/src/McSolver/Libraries/Simulator/Logic/Routines/StatisticsRoutines.c +++ b/src/McSolver/Libraries/Simulator/Logic/Routines/StatisticsRoutines.c @@ -149,6 +149,15 @@ Vector3_t CalculateMobilityVector(SCONTEXT_PARAMETER, const Vector3_t *displacem return result; } +Vector3_t CalculateConductivityVector(SCONTEXT_PARAMETER, const Vector3_t *mobility, double charge, double particleDensity) +{ + Vector3_t result; + result.A = CalculateTotalConductivity(mobility->A, charge, particleDensity); + result.B = CalculateTotalConductivity(mobility->B, charge, particleDensity); + result.C = CalculateTotalConductivity(mobility->C, charge, particleDensity); + return result; +} + double CalculateParticleDensity(SCONTEXT_PARAMETER, byte_t particleId) { let volume = CalculateSuperCellVolume(simContext); @@ -186,6 +195,15 @@ static inline double CalculateAverageMigrationRate(double simulatedTime, int64_t Vector3_t CalculateNernstEinsteinConductivity(SCONTEXT_PARAMETER, const Vector3_t *diffusionVector, const byte_t particleId, const double particleDensity) { let charge = getDbStructureModelMetaData(simContext)->ParticleCharges[particleId] * NATCONST_ELMCHARGE; + let normValues = CalculateNormalizedNernstEinsteinConductivity(simContext, diffusionVector, particleDensity); + + let result = ScalarMultiplyVector3(&normValues, charge * charge); + return result; +} + +Vector3_t CalculateNormalizedNernstEinsteinConductivity(SCONTEXT_PARAMETER, const Vector3_t *diffusionVector, double particleDensity) +{ + let charge = 1.0; let tempFactor = getDbModelJobInfo(simContext)->Temperature * NATCONST_BLOTZMANN * NATCONST_ELMCHARGE; let factor = charge * charge * particleDensity / tempFactor; @@ -193,6 +211,15 @@ Vector3_t CalculateNernstEinsteinConductivity(SCONTEXT_PARAMETER, const Vector3_ return result; } +Vector3_t CalculateDiffusionCoefficientsSigma(SCONTEXT_PARAMETER, const Vector3_t *conductivities, byte_t particleId, double particleDensity) +{ + let charge = getDbStructureModelMetaData(simContext)->ParticleCharges[particleId] * NATCONST_ELMCHARGE; + let tempFactor = getDbModelJobInfo(simContext)->Temperature * NATCONST_BLOTZMANN * NATCONST_ELMCHARGE; + let factor = tempFactor / (particleDensity * charge * charge); + let result = ScalarMultiplyVector3(conductivities, factor); + return result; +} + void PopulateMobilityData(SCONTEXT_PARAMETER, ParticleMobilityData_t *restrict data) { let meta = getDbStructureModelMetaData(simContext); @@ -214,13 +241,15 @@ void PopulateMobilityData(SCONTEXT_PARAMETER, ParticleMobilityData_t *restrict d vector3ScalarOp(data->MeanMoveR2, count, /=); data->MigrationRate = CalculateAverageMigrationRate(simulatedTime, data->ParticleStatistics->CounterCollection->McsCount, count); - data->DiffusionCoefficient = CalculateDiffusionCoefficient(&data->MeanMoveR2, simulatedTime); + data->DiffusionCoefficientVector = CalculateDiffusionCoefficient(&data->MeanMoveR2, simulatedTime); data->MobilityVector = CalculateMobilityVector(simContext, &data->MeanMoveR1, &meta->NormElectricFieldVector); data->TotalMobility = CalculateFieldProjectedMobility(simContext, &data->MeanMoveR1, &meta->NormElectricFieldVector); - data->NernstEinsteinConductivity = CalculateNernstEinsteinConductivity(simContext, &data->DiffusionCoefficient, id, density); + data->ConductivityVector = CalculateConductivityVector(simContext, &data->MobilityVector, meta->ParticleCharges[id], density); data->TotalConductivity = CalculateTotalConductivity(data->TotalMobility, meta->ParticleCharges[id], density); - data->TotalConductivityPerCharge = CalculateTotalConductivity(data->TotalMobility, 1, density); + data->NormalizedConductivityVector = CalculateConductivityVector(simContext, &data->MobilityVector, 1.0, density); + + data->DSigmaVector = CalculateDiffusionCoefficientsSigma(simContext, &data->ConductivityVector, id, density); } void PopulateParticleStatistics(SCONTEXT_PARAMETER, ParticleStatistics_t *restrict statistics) diff --git a/src/McSolver/Libraries/Simulator/Logic/Routines/StatisticsRoutines.h b/src/McSolver/Libraries/Simulator/Logic/Routines/StatisticsRoutines.h index f2f8424c..e04c0515 100644 --- a/src/McSolver/Libraries/Simulator/Logic/Routines/StatisticsRoutines.h +++ b/src/McSolver/Libraries/Simulator/Logic/Routines/StatisticsRoutines.h @@ -47,20 +47,23 @@ typedef struct ParticleMobilityData // The total mobility in [m^2/(V s)] double TotalMobility; - // The total conductivity in [S/m] + // The total conductivity in [S/m] in field direction double TotalConductivity; - // The normalized conductivity per positive charge number [S/(m*z)] - double TotalConductivityPerCharge; - // The actual migration rate in [Hz] double MigrationRate; // The total diffusion coefficient components in [m^2/s] in (x,y,z) direction - Vector3_t DiffusionCoefficient; + Vector3_t DiffusionCoefficientVector; + + // The sigma diffusion coefficient components in [m^2/s] in (x,y,z) direction as calculated from the conductivity + Vector3_t DSigmaVector; - // The conductivity components in [S/m] as calculated using the Stokes-Einstein relation (Fluctuation dissipation theorem) - Vector3_t NernstEinsteinConductivity; + // The conductivity vector in [m^2/(V s)] + Vector3_t ConductivityVector; + + // The conductivity vector in [m^2/(V s)] for a pseudo charge of +1 + Vector3_t NormalizedConductivityVector; // The mobility vector in [m^2/(V s)] Vector3_t MobilityVector; @@ -186,10 +189,20 @@ Vector3_t GetGlobalTrackerEnsembleShift(SCONTEXT_PARAMETER, JumpCollection_t *ju // Calculates a mobility vector in [m^2/(V s)] using the provided mean displacement and normalized field vector Vector3_t CalculateMobilityVector(SCONTEXT_PARAMETER, const Vector3_t *displacement, const Vector3_t *normFieldVector); +// Calculates a conductivity vector in [m^2/(V s)] using the provided mean displacement and particle density and charge z +Vector3_t CalculateConductivityVector(SCONTEXT_PARAMETER, const Vector3_t *mobility, double charge, double particleDensity); + // Calculates the conductivity components (x,y,z) for the passed particle id in [S/m] from diffusion components and particle density -// using the fluctuation dissipation theorem +// using the nernst einstein relation Vector3_t CalculateNernstEinsteinConductivity(SCONTEXT_PARAMETER, const Vector3_t *diffusionVector, byte_t particleId, double particleDensity); +// Calculates the conductivity components (x,y,z) for the passed particle id in [S/m] from diffusion components and particle density +// using the nernst einstein relation (normalized to a charge of +1) +Vector3_t CalculateNormalizedNernstEinsteinConductivity(SCONTEXT_PARAMETER, const Vector3_t *diffusionVector, double particleDensity); + +// Calculates the diffusion coefficients [m^2 / s] as defined by the Nernst-Einstein relation +Vector3_t CalculateDiffusionCoefficientsSigma(SCONTEXT_PARAMETER, const Vector3_t *conductivities, byte_t particleId, double particleDensity); + // Calculates the total mobility in [m^2/(V s)] using the provided mean displacement and normalized field vector double CalculateFieldProjectedMobility(SCONTEXT_PARAMETER, const Vector3_t *displacement, const Vector3_t *normFieldVector); diff --git a/src/ModelBuilder/Mocassin.UI.GUI/Controls/DxVisualizer/ModelViewer/DxModelSceneViewModel.cs b/src/ModelBuilder/Mocassin.UI.GUI/Controls/DxVisualizer/ModelViewer/DxModelSceneViewModel.cs index fdcd50bb..d50e0d27 100644 --- a/src/ModelBuilder/Mocassin.UI.GUI/Controls/DxVisualizer/ModelViewer/DxModelSceneViewModel.cs +++ b/src/ModelBuilder/Mocassin.UI.GUI/Controls/DxVisualizer/ModelViewer/DxModelSceneViewModel.cs @@ -557,8 +557,8 @@ private void AddCellFrameSceneNode(DxSceneBuilder sceneBuilder, in FractionalBox var matrix = VectorTransformer.FractionalSystem.ToCartesianMatrix.ToDxMatrix(); var lineBuilder = new LineBuilder(); - var (aMax, bMax, cMax) = (renderBox.End.A.FloorToInt(), renderBox.End.B.FloorToInt(), renderBox.End.C.FloorToInt()); - var (aMin, bMin, cMin) = (renderBox.Start.A.CeilToInt(), renderBox.Start.B.CeilToInt(), renderBox.Start.C.CeilToInt()); + var (aMax, bMax, cMax) = (renderBox.End.A.CeilToInt(), renderBox.End.B.CeilToInt(), renderBox.End.C.CeilToInt()); + var (aMin, bMin, cMin) = (renderBox.Start.A.FloorToInt(), renderBox.Start.B.FloorToInt(), renderBox.Start.C.FloorToInt()); for (var a = aMin; a <= aMax; a++) { for (var b = bMin; b <= bMax; b++) diff --git a/src/ModelBuilder/Mocassin.UI.Xml/StructureModel/CellReferencePositionData.cs b/src/ModelBuilder/Mocassin.UI.Xml/StructureModel/CellReferencePositionData.cs index 66c54244..0ff03802 100644 --- a/src/ModelBuilder/Mocassin.UI.Xml/StructureModel/CellReferencePositionData.cs +++ b/src/ModelBuilder/Mocassin.UI.Xml/StructureModel/CellReferencePositionData.cs @@ -77,7 +77,7 @@ protected override ModelObject GetModelObjectInternal() { OccupationSet = new ParticleSet {Key = Occupation.Key}, Stability = Stability, - Vector = new Fractional3D(Math.Abs(A), Math.Abs(B), Math.Abs(C)) + Vector = new Fractional3D(A, B, C).TrimToUnitCell(1e-13) }; return obj; }