Autor |
Nachricht |
Lorz
Anmeldungsdatum: 18.09.2014 Beiträge: 74
|
|
|
TomS Moderator
Anmeldungsdatum: 20.03.2009 Beiträge: 18090
|
TomS Verfasst am: 10. Feb 2024 14:38 Titel: |
|
|
Man muss ein gekoppeltes lineares Gleichungssystem ansetzen, das alle Zerfallskanäle berücksichtigt.
Wenn man immer nur einen Zerfallskanal je Nuklid hat, wird das ziemlich einfach, da man eine DGL nach der anderen lösen kann.
_________________ Niels Bohr brainwashed a whole generation of theorists into thinking that the job (interpreting quantum theory) was done 50 years ago. |
|
|
Lorz
Anmeldungsdatum: 18.09.2014 Beiträge: 74
|
Lorz Verfasst am: 10. Feb 2024 15:06 Titel: |
|
|
@TomS
Danke für die Antwort. Dachte ich mir schon, dass ich das nicht so einfach ansetzen kann. Bzw, das war mein erster Gedanke. Danach habe ich gedacht - ist ja doch nicht kompliziert. Aber der erste Gedanke war wohl schon richtig.
Konkret zu Ergebnissen:
Ist die Anzahl noch korrekt? Denn diese Art von Kurve finde ich auch in der Literatur.
Von der Berechnung beliebig vieler Generationen sehe ich dann gerne ab. Aber kann mir jemand die Funktion noch für die dritte Generation nennen?
Ein bisschen versuche ich mitzuhelfen:
Also der zukünftige Bestand (zur Zeit ) der Generation 3 ergibt sich aus dem Zuwachs durch den momentanen Abbau von Generation 2 (gleichbedeutend mit Umwandlung in die Generation 3) innerhalb zur Zeit und dem Abbau des schon vorhandenen Bestands von :
Soweit korrekt?
Da der Anfangswert für die erste Generation als gegeben vorausgesetzt ist, könnte ich meinem Verständnis nach, die 2. Generation berechnen, denn und sind ja bekannt. Und iterativ sollten dann auch alle anderen Generationen berechenbar sein.
Aber ich weiß halt nicht mehr, wie das geht. Also ich kann mich an so ein paar Verfahren (Picard-Lindelöff oder so) erinnern aus der Analysis. Theo-Physik ist leider fast komplett vergessen bei mir.
Also wenn mir jemand einfach so die Lösung für Generation 3 sagen kann, bin ich happy :-)
|
|
|
3.Generation Gast
|
3.Generation Verfasst am: 10. Feb 2024 16:50 Titel: |
|
|
Lorz hat Folgendes geschrieben: |
Also wenn mir jemand einfach so die Lösung für Generation 3 sagen kann, bin ich happy :-) |
Das funktioniert nur wenn
Bei wird die Sache wesentlich komplizierter
|
|
|
Lorz
Anmeldungsdatum: 18.09.2014 Beiträge: 74
|
Lorz Verfasst am: 10. Feb 2024 16:55 Titel: |
|
|
@ Gast Danke!
Ich verstehe allerdings gerade nicht, wieso bereits bei
ein auftaucht.
Auch bin ich skeptisch, was die Dimension/Einheiten der Summanden in und betrifft. Also die Zerfallskonstante ist ja nicht einfach dimensionslos (hat doch ) . Oder wie verhält sich dies?
|
|
|
3.Generation Gast
|
3.Generation Verfasst am: 10. Feb 2024 18:01 Titel: |
|
|
Lorz hat Folgendes geschrieben: |
Ich verstehe allerdings gerade nicht, wieso bereits bei
ein auftaucht.
|
Man löst diese Aufgabe am besten mit einer Eigenwertrechnung
und der Eigenvektor von ist
Es gilt zu jedem Zeitpunkt
Lorz hat Folgendes geschrieben: |
Auch bin ich skeptisch, was die Dimension/Einheiten der Summanden in und betrifft. Also die Zerfallskonstante ist ja nicht einfach dimensionslos (hat doch ) . Oder wie verhält sich dies? |
Gut erkannt
Man könnte die Einheit ausgleichen indem man s bei C1 dazuschreibt
|
|
|
Lorz
Anmeldungsdatum: 18.09.2014 Beiträge: 74
|
Lorz Verfasst am: 10. Feb 2024 18:16 Titel: |
|
|
Ich probiere jetzt mal selbst.
Also die Gleichung
setze ich jetzt erst mal für die zweite Generation an und schaue, ob das was vernünftiges bei rauskommt.
Umformungen
...ähhm, ja...
|
|
|
Lorz
Anmeldungsdatum: 18.09.2014 Beiträge: 74
|
Lorz Verfasst am: 10. Feb 2024 18:22 Titel: |
|
|
3.Generation hat Folgendes geschrieben: | Lorz hat Folgendes geschrieben: |
Ich verstehe allerdings gerade nicht, wieso bereits bei
ein auftaucht.
|
Man löst diese Aufgabe am besten mit einer Eigenwertrechnung
und der Eigenvektor von ist
Es gilt zu jedem Zeitpunkt
Lorz hat Folgendes geschrieben: |
Auch bin ich skeptisch, was die Dimension/Einheiten der Summanden in und betrifft. Also die Zerfallskonstante ist ja nicht einfach dimensionslos (hat doch ) . Oder wie verhält sich dies? |
Gut erkannt
Man könnte die Einheit ausgleichen indem man s bei C1 dazuschreibt |
OK, verstehe, man hat hier viel mehr Flexibilität. Ich kann mit einer beliebigen Anfangssituation starten.
Bei meinem vereinfachten Beispiel (alle Mutternuklide warten mit dem Zerfallen bis , dann gilt
Weiter müsste nach der obersten Gleichung gelten
Soll denn die Zerfallskonstante sein? Und entsprechend ?
Dann wäre
Mit erhalten
Und somit für
Ich prüfe mal für mich die Plausibität von :
Also für hebt sich alles weg, , das passt.
Nun nehme ich mal den Fall, dass winzig klein ist, dann dominierend gegenüber in der Differenz und kürzt sich im Bruch weg. Wir erhalten
bzw
Auch das macht Sinn, die Mutternuklide wandeln sich quasi sofort um, und
und nur noch
bestimmt den Verlauf.
Was ich nicht verstehe - warum "explodiert" die Anzahl, wenn ?
Und warum kommt in der Lösung für kein vor?
|
|
|
3.Generation Gast
|
3.Generation Verfasst am: 10. Feb 2024 19:44 Titel: |
|
|
Lorz hat Folgendes geschrieben: |
Was ich nicht verstehe - warum "explodiert" die Anzahl, wenn ?
|
Da explodiert nichts,weil die Gesamtmenge hier immer N1(0) beträgt
Man hat hier eine Grenzwertrechnung
So etwa
Lorz hat Folgendes geschrieben: |
Und warum kommt in der Lösung für kein vor? |
Streng genommen kommt vor
Wenn N3 nicht weiter zerfällt ist und damit
|
|
|
TomS Moderator
Anmeldungsdatum: 20.03.2009 Beiträge: 18090
|
TomS Verfasst am: 10. Feb 2024 20:03 Titel: |
|
|
Man setzt:
Im der k-ten Zeile steht die Änderung des k-ten Nuklids, in der Klammer mit dem Minuszeichen die Zerfälle in die folgenden Tochternuklide k+1, k+2 … , in der Klammer mit dem Pluszeichen die Zerfälle aus den Mutternukliden k-1, k-2 …
Das kann man als Gleichung mit einer Zerfallsmatrix Lambda und einem Vektor N auffassen:
Auf den Diagonalelementen kk stehen die Summen der Zerfallskonstanten in die Tochternuklide. In einem Dreieck die einzelnen Matrixelemente km mit m>k mit den Zerfallskonstanten von k nach m. Das zweite Dreieck ist Null, keine Zerfälle von m nach k für m>k.
In einfachen Fällen d.h. für wenige beteiligte Nuklide bzw. viele Nullen man das mittels Diagonalisierung von Lambda analytisch lösen. Ansonsten füttert man das in einen Computer.
_________________ Niels Bohr brainwashed a whole generation of theorists into thinking that the job (interpreting quantum theory) was done 50 years ago. |
|
|
Lorz
Anmeldungsdatum: 18.09.2014 Beiträge: 74
|
Lorz Verfasst am: 10. Feb 2024 20:39 Titel: |
|
|
3.Generation hat Folgendes geschrieben: | Lorz hat Folgendes geschrieben: |
Was ich nicht verstehe - warum "explodiert" die Anzahl, wenn ?
|
Da explodiert nichts,weil die Gesamtmenge hier immer N1(0) beträgt
Man hat hier eine Grenzwertrechnung
So etwa
|
Ah, ich habe nicht Summe bzw Differenz aus beiden Brüche betrachtet.
Einzeln gesehen geht der Bruch ja schon gegen unendlich. Ist daher hier nicht das Grenzwertproblem, dass
entsteht?
3.Generation hat Folgendes geschrieben: |
Lorz hat Folgendes geschrieben: |
Und warum kommt in der Lösung für kein vor? |
Streng genommen kommt vor
Wenn N3 nicht weiter zerfällt ist und damit |
Kann man denn noch angenehm erweitern auf eine Halbwertszeit 3. Generation kleiner als unendlich?
Zuletzt bearbeitet von Lorz am 10. Feb 2024 20:44, insgesamt einmal bearbeitet |
|
|
Lorz
Anmeldungsdatum: 18.09.2014 Beiträge: 74
|
Lorz Verfasst am: 10. Feb 2024 20:44 Titel: |
|
|
TomS hat Folgendes geschrieben: | Man setzt:
Im der k-ten Zeile steht die Änderung des k-ten Nuklids, in der Klammer mit dem Minuszeichen die Zerfälle in die folgenden Tochternuklide k+1, k+2 … , in der Klammer mit dem Pluszeichen die Zerfälle aus den Mutternukliden k-1, k-2 …
Das kann man als Gleichung mit einer Zerfallsmatrix Lambda und einem Vektor N auffassen:
Auf den Diagonalelementen kk stehen die Summen der Zerfallskonstanten in die Tochternuklide. In einem Dreieck die einzelnen Matrixelemente km mit m>k mit den Zerfallskonstanten von k nach m. Das zweite Dreieck ist Null, keine Zerfälle von m nach k für m>k.
In einfachen Fällen d.h. für wenige beteiligte Nuklide bzw. viele Nullen man das mittels Diagonalisierung von Lambda analytisch lösen. Ansonsten füttert man das in einen Computer. |
OK, damit könnte ich ja tatsächlich beliebig lange Zerfallsreihen behandeln.
Aber ich verstehen nicht den Doppelindex, an der Zerfallskonstanten. Nach meiner Idee habe ich für jedes der k Generationen an Nukliden eine Zerfallskonstante . Ich kann diese vielen nicht zuordnen.
|
|
|
TomS Moderator
Anmeldungsdatum: 20.03.2009 Beiträge: 18090
|
TomS Verfasst am: 10. Feb 2024 20:50 Titel: |
|
|
Wenn’s dich interessiert, hier der Einzeiler in Python:
Code: | qty_sol_space[:,i] = np.dot(la.expm(t_space[i] * dec_mat), qty_nuclides_0) |
Das löst das DGL-System für beliebige Zeiten bei gegebener Anfangsbedingung.
Außenrum noch 50 Zeilen Code zum Befüllen von Lambda und für die Graphik.
Ich hätte das für die Zerfallsreihe von Radium fertig und würde es dir nach dem Urlaub schicken.
_________________ Niels Bohr brainwashed a whole generation of theorists into thinking that the job (interpreting quantum theory) was done 50 years ago. |
|
|
TomS Moderator
Anmeldungsdatum: 20.03.2009 Beiträge: 18090
|
TomS Verfasst am: 10. Feb 2024 20:54 Titel: |
|
|
Lorz hat Folgendes geschrieben: | Aber ich verstehen nicht den Doppelindex, an der Zerfallskonstanten. |
Im Allgemeinen kann ein Nuklid k in verschiedene Tochternuklide zerfallen, also in m1, m2 … – siehe Abbildung.
Ist das nicht der Fall, bzw. setzt man näherungsweise an, dass je Nuklid ein Zerfall dominiert, d.h. dass eine Zerfallsreihe vorliegt, dann hast du immer nur zwei Einträge für einen Zerfall:
km mit k<m: Nuklid k wird weniger
mk mit k<m: Nuklid m wird um den selben Betrag mehr
Beide Einträge in der Matrix haben unterschiedliches Vorzeichen, d.h. letztlich kommt je Nuklid k nur eine Zahl vor, die Zerfallskonstante k. Und diese steht genau zweimal da, dämlich bei
km mit m=k+1: Nuklid k wird weniger
mk mit m=k+1: Nuklid m=k+1 wird um den selben Betrag mehr
Beschreibung: |
|
Download |
Dateiname: |
IMG_3805.png |
Dateigröße: |
144.94 KB |
Heruntergeladen: |
22 mal |
_________________ Niels Bohr brainwashed a whole generation of theorists into thinking that the job (interpreting quantum theory) was done 50 years ago. |
|
|
TomS Moderator
Anmeldungsdatum: 20.03.2009 Beiträge: 18090
|
|
|
Lorz
Anmeldungsdatum: 18.09.2014 Beiträge: 74
|
Lorz Verfasst am: 11. Feb 2024 10:02 Titel: |
|
|
TomS hat Folgendes geschrieben: | Lorz hat Folgendes geschrieben: | Aber ich verstehen nicht den Doppelindex, an der Zerfallskonstanten. |
Im Allgemeinen kann ein Nuklid k in verschiedene Tochternuklide zerfallen, also in m1, m2 … – siehe Abbildung.
Ist das nicht der Fall, bzw. setzt man näherungsweise an, dass je Nuklid ein Zerfall dominiert, d.h. dass eine Zerfallsreihe vorliegt, dann hast du immer nur zwei Einträge für einen Zerfall:
km mit k<m: Nuklid k wird weniger
mk mit k<m: Nuklid m wird um den selben Betrag mehr
Beide Einträge in der Matrix haben unterschiedliches Vorzeichen, d.h. letztlich kommt je Nuklid k nur eine Zahl vor, die Zerfallskonstante k. Und diese steht genau zweimal da, dämlich bei
km mit m=k+1: Nuklid k wird weniger
mk mit m=k+1: Nuklid m=k+1 wird um den selben Betrag mehr |
Ah, OK, verstehe, mehrere Zerfallsmöglichkeiten pro Nuklid, mehrere Zerfallskonstanten.
Danke für das Bild, ja, klar, Zerfallsreihen mit mehr als einem Weg kenne ich.
Aber nu mach mal ruhig Deinen Urlaub! Ich probiere mich mal an Deine dargebotenen Möglichkeiten, sollte ja recht knapp werden meine Matrix, wenn ich erstmal nur die bis ink. der 3. Generation gehe und dabei nur eine Zerfallskonstante pro Radionuklid annehme.
Ansonsten interessiert mich schon die weitere Entwicklung der Gesamtaktivität. Also mein ursprüngliches Ziel war es, zu verstehen, wie es zu einem säkularen Gleichgewicht kommt, bei dem das konstante Plateau recht lang ist.
|
|
|
TomS Moderator
Anmeldungsdatum: 20.03.2009 Beiträge: 18090
|
TomS Verfasst am: 11. Feb 2024 13:40 Titel: |
|
|
Die Gesamtaktivität für alpha-und beta-Strahlung purzelt auch heraus; gamma-Strahlung nicht, das müsste man separat lösen.
_________________ Niels Bohr brainwashed a whole generation of theorists into thinking that the job (interpreting quantum theory) was done 50 years ago. |
|
|
TomS Moderator
Anmeldungsdatum: 20.03.2009 Beiträge: 18090
|
TomS Verfasst am: 19. Feb 2024 09:25 Titel: |
|
|
Ich habe jetzt mal die vollständige Zerfallskette
mit dem Differentialgleichungssystem
für den Vektor N nach zwei Methoden berechnet, 1. numerische Lösung, 2. Matrix-Exponential
Beide Methoden liefern bis auf numerische Fehler < 1e-3 übereinstimmende Ergebnisse.
Anbei zunächst die Anzahl der jeweiligen Nuklide N, die Aktivitäten
sowie die freiwerdenden Energien
Als Anfangsbedingung setzte ich 1 mol Ra226 an.
Hier der Python Code
Code: | #
""" decay chain;
example: Radium-226
attetion: works for linear decay chains only
Author: Tom
Date: 2024-02-15
"""
from enum import Enum, auto
import numpy as np
import scipy.constants as sp_sc
import scipy.linalg as sp_la
import scipy.integrate as sp_in
import matplotlib.pyplot as plt
plt.rc('font', size=8)
import mplcyberpunk
plt.style.use("cyberpunk")
# the solution ------------------------------
# using the matrix exponential: calculate the quantity, activity and energy for all nuclides as a function of time
def calc_dec_chain(Nuclides, dec_coeffs, dec_energies, qty_nuclides_0, t_start, t_end, t_space):
def init_dec_mat(Nuclides, dec_coeffs):
""" restricted to linear decay chains """
dec_mat = np.zeros((len(Nuclides),len(Nuclides)))
for m in range (len(Nuclides)):
for n in range (len(Nuclides)):
if (m == n):
dec_mat[m,n] = -dec_coeffs[n]
if (n > 0):
dec_mat[m,n-1] = dec_coeffs[n-1]
return dec_mat
dec_mat = init_dec_mat(Nuclides, dec_coeffs)
qty_sol_space = np.zeros((len(Nuclides),t_nsamples))
act_sol_space = np.zeros((len(Nuclides),t_nsamples))
eny_sol_space = np.zeros((len(Nuclides),t_nsamples))
for i in range(t_nsamples):
qty_sol_space[:,i] = np.dot(sp_la.expm(t_space[i] * dec_mat), qty_nuclides_0)
act_sol_space[:,i] = dec_coeffs * qty_sol_space[:,i]
eny_sol_space[:,i] = dec_energies * dec_coeffs * qty_sol_space[:,i]
return qty_sol_space, act_sol_space, eny_sol_space
# using a 1st order system of differential equations: calculate the quantity, activity and energy for all nuclides as a function of time
def solve_dec_chain(Nuclides, dec_coeffs, dec_energies, qty_nuclides_0, t_start, t_end, t_space):
def def_decay_chain(t, q, c):
""" restricted to linear decay chains """
dot_q = np.zeros_like(q)
dot_q[0] = -c[0] * q[0]
dot_q[-1] = c[-2] * q[-2]
for n in range(1, len(q)-1):
dot_q[n] = c[n-1] * q[n-1] - c[n] * q[n]
return dot_q
t_span = (t_start, t_end)
sol = sp_in.solve_ivp(def_decay_chain, t_span, qty_nuclides_0, args=(dec_coeffs,), t_eval=t_space, method='BDF')
return sol.y, None, None
# extract alpha and beta channels, calculate sums for alpha and beta
def sum_alpha_beta_channels(sol_space):
alpha = np.zeros(t_nsamples)
beta = np.zeros(t_nsamples)
for i in range(t_nsamples):
alpha[i] += (sol_space[Nuclides.Ra_226.value, i] +
sol_space[Nuclides.Rn_222.value, i] +
sol_space[Nuclides.Po_218.value, i] +
sol_space[Nuclides.Po_214.value, i] +
sol_space[Nuclides.Po_210.value, i])
beta[i] += (sol_space[Nuclides.Pb_214.value, i] +
sol_space[ Nuclides.Bi_214.value, i] +
sol_space[ Nuclides.Pb_210.value, i] +
sol_space[ Nuclides.Bi_210.value, i])
return alpha, beta
# plot functions ------------------------------
def subplot_quantities(t_space_years, qty_sol_space):
for n in Nuclides:
plt.plot(t_space_years, qty_sol_space[n.value], lw=0.75, ls='-', label=n.name)
return
def subplot_activities(t_space_years, qty_sol_space):
for n in Nuclides:
if (n.value == len(Nuclides) - 1):
break
plt.plot(t_space_years, act_sol_space[n.value], lw=0.75, ls='-', label=n.name)
return
def subplot_energies(t_space_years, qty_sol_space):
for n in Nuclides:
if (n.value == len(Nuclides) - 1):
break
plt.plot(t_space_years, eny_sol_space[n.value], lw=0.75, ls='-', label=n.name)
return
def subplot_channels(t_space_years, radiation):
plt.plot(t_space_years, radiation[0], lw=0.75, ls='-', label='alpha')
plt.plot(t_space_years, radiation[1], lw=0.75, ls='-', label='beta')
return
def plot_frame(t_space, qty_sol_space, x_label, y_label, title, subplot_func):
plt.figure()
plt.xscale('log')
plt.yscale('log')
plt.xlabel(x_label)
plt.ylabel(y_label)
plt.title(title)
subplot_func(t_space/sp_sc.year, qty_sol_space)
plt.legend()
plt.show()
return
# the nuclides in the Ra-226 decay chain, their decay half-times [seconds], their decay energies [MeV], and their quantity ------------------------------
Nuclides = Enum('Nuclides', ['Ra_226', 'Rn_222', 'Po_218', 'Pb_214', 'Bi_214', 'Po_214', 'Pb_210', 'Bi_210', 'Po_210', 'Pb_206'], start=0)
dec_times = np.array([1600 * sp_sc.year, 3.82 * sp_sc.day, 3.098 * sp_sc.minute, 27.06 * sp_sc.minute, 19.9 * sp_sc.minute, 0.0001643, 22.20 * sp_sc.year, 5.012 * sp_sc.day, 138.376 * sp_sc.day, np.inf])
dec_coeffs = np.log(2) / dec_times
dec_energies = np.array([4.7, 5.5, 6.0, 0.3, 1.0, 7.7, 0.006, 0.3, 5.3, 0.0])
qty_nuclides_0 = np.zeros_like(dec_times)
qty_nuclides_0[0] = sp_sc.Avogadro
# driver code ------------------------------
# set time span
t_start = 1.0 * sp_sc.hour
t_end = 10000.0 * sp_sc.year
t_nsamples = 1000
t_space = np.logspace(np.log10(t_start), np.log10(t_end), t_nsamples)
# solve equations
qty_sol_space, act_sol_space, eny_sol_space = calc_dec_chain(Nuclides, dec_coeffs, dec_energies, qty_nuclides_0, t_start, t_end, t_space)
qty_sol_space_2, _, _ = solve_dec_chain(Nuclides, dec_coeffs, dec_energies, qty_nuclides_0, t_start, t_end, t_space)
# sum for alpha and beta channels
act_alpha, act_beta = sum_alpha_beta_channels(act_sol_space)
eny_alpha, eny_beta = sum_alpha_beta_channels(eny_sol_space)
# plot
plot_frame(t_space, qty_sol_space, 'time [years]', 'quantities [N]','quantity', subplot_quantities)
plot_frame(t_space, qty_sol_space_2, 'time [years]', 'quantities [N]','quantity', subplot_quantities)
plot_frame(t_space, qty_sol_space_2/qty_sol_space-1.0, 'time [years]', 'relative error','relative error', subplot_quantities)
plot_frame(t_space, qty_sol_space, 'time [years]', 'activities [Bq]','activity', subplot_activities)
plot_frame(t_space, qty_sol_space, 'time [years]', 'energies [MeV]','energy', subplot_energies)
plot_frame(t_space, (act_alpha, act_beta), 'time [years]', 'activities [Bq]','activity', subplot_channels)
plot_frame(t_space, (eny_alpha, eny_beta), 'time [years]', 'energies [MeV]','energy', subplot_channels) |
Beschreibung: |
|
Dateigröße: |
43.26 KB |
Angeschaut: |
556 mal |
|
Beschreibung: |
|
Dateigröße: |
33.14 KB |
Angeschaut: |
556 mal |
|
Beschreibung: |
|
Dateigröße: |
48.35 KB |
Angeschaut: |
556 mal |
|
_________________ Niels Bohr brainwashed a whole generation of theorists into thinking that the job (interpreting quantum theory) was done 50 years ago.
Zuletzt bearbeitet von TomS am 20. Feb 2024 22:32, insgesamt 3-mal bearbeitet |
|
|
TomS Moderator
Anmeldungsdatum: 20.03.2009 Beiträge: 18090
|
|
|
|
|