Scientia et Technica Año XXVIII, Vol. 29, No. 04, octubre-diciembre de 2024. Universidad Tecnológica de Pereira. ISSN 0122-1701 y ISSN-e: 2344-7214 167

Abstract— Nuclear reactivity and neutron population

density are critical parameters for quantifying and

describing the behavior of a nuclear reactor. Accurate

knowledge of reactivity is essential for the safe control and

operation of a reactor. Various approaches and methods

for calculating reactivity are reported in the literature. In

this paper, we propose using the Euler-Maclaurin formula

to numerically solve reactivity through the inverse point

kinetics equation, a key model in the design of digital

reactivity meters. By approximating an integral as a sum

in the expansion of a continuous integral into a discrete

version, and considering only the first three Bernoulli

numbers while reducing higher-order derivatives, we

derive a model to calculate nuclear reactivity as a function

of neutron population density and physical constants in

thermal reactors using uranium as nuclear fuel.

Index Terms— Bernoulli numbers; Euler-Maclaurin formula;

Neutron density; Nuclear reactivity; Nuclear reactor; Numerical

simulation.

Resumen—La reactividad nuclear y la densidad de la población

de neutrones son los parámetros más significativos al momento

de cuantificar y describir físicamente el comportamiento de un

reactor nuclear. Conocer el valor de la reactividad con suficiente

precisión permite realizar el control y la operación segura de un

reactor nuclear. Existen varios enfoques y métodos para calcular

la reactividad reportados en la literatura. Sin embargo, en este

artículo se propone usar la fórmula de Euler-Maclaurin para la

solución numérica de la reactividad mediante la ecuación inversa

de la cinética puntual que ha sido un modelo clave en el diseño de

medidores de reactividad digital. Considerando la aproximación

de una integral como una suma, en el desarrollo de la expansión

de la integral continua en una versión discreta y tomando en

cuenta únicamente los tres primeros números de Bernoulli y una

Este manuscrito fue enviado el día 11 agosto, 2023, aceptado el 10 de

septiembre de 2024 y publicado el 20 de diciembre del 2024

Este artículo fue desarrollado con el apoyo financiero de la Universidad

Surcolombiana.

Daniel Suescún Díaz es un investigador del grupo de Física aplicada,

FIASUR de la Universidad surcolombiana, Neiva, Colombia. (e-mail:

daniel.suescun@usco.edu.co).

Geraldyne Ule Duque es una investigadora del grupo de Física aplicada,

FIASUR de la Universidad surcolombiana, Neiva, Colombia. (e-mail:

geraldinu525@gmail.com).

Diego Peña Lara es un investigador del grupo de Transiciones de Fase y

Materiales Funcionales de la Universidad del Valle, Cali, Colombia. (e-mail:

diego.pena@correounivalle.edu.co).

aproximación que reduce las derivadas de mayor orden, es

posible deducir un modelo que permite el cálculo de la

reactividad nuclear en función de la densidad de la población de

neutrones y de las constantes física en reactores térmicos que

usan uranio como combustible nuclear.

Palabras claves— Números de Bernoulli; Formula Euler-

Maclaurin; Densidad de neutrones; reactividad nuclear; Reactor

nuclear; Simulación numérica.

I. INTRODUCCIÓN

no de los retos de la actualidad es producir electricidad para el

desarrollo económico de los pueblos sin aumentar el

calentamiento global del planeta, un camino para lograrlo es

mediante la energía nuclear que se obtiene a partir de los

núcleos atómicos en un proceso conocido como fisión nuclear.

En esta reacción nuclear se liberan alrededor de 200 MeV con

una producción de varios neutrones, los cuales van a producir

nuevas fisiones. Para una reacción en cadena controlada es

necesario un reactor nuclear [1]. Para dicho control, es

necesario conocer con mayor precisión un parámetro llamado

reactividad obtenida a partir de una forma conocida de la

densidad de neutrones.

Diferentes métodos experimentales y numéricos se utilizan

para estimar y calcular la reactividad, este valor se calcula con

la ecuación inversa de la cinética puntual, la cual es una

ecuación integro diferencial que tiene asociada a la parte

diferencial el periodo del reactor y la parte integral contiene

todos los valores de la densidad de la población de neutrones.

A lo largo de la historia de los reactores nucleares se han

llevado a cabo diferentes estudios para calcular la reactividad,

un trabajo pionero usa diferencias finitas con el filtro pasa-

bajo [2]. Otro trabajo usa un método estadístico para

determinar la fuente de neutrones [3]. Posteriormente, se

presenta un método que puede calcular la reactividad en

tiempo real con una aproximación lineal [4]. Análogamente,

se propone una medición de la reactividad en tiempo real [5].

El método conocido como memorial índex se presenta para el

cálculo de la reactividad [6]. Así mismo, se publicó un nuevo

método usando los polinomios de Lagrange de tercero y

quinto orden de precisión [7]. La serie de Euler-Maclaurin es

usada para el cálculo de la reactividad con una buena

aproximación [8]. Otro método para el cálculo de la

reactividad emplea un pronóstico-corrección de alta precisión

Hamming [9]. En otros métodos se hace uso de la

Solución Numérica de la Reactividad en

Reactores Nucleares

Numerical Solution of Reactivity for Nuclear Reactors

D. Suescún-Díaz ; G. Ule-Duque ; D. Peña-Lara

DOI: https://doi.org/10.22517/23447214.25438

Scientific and technological research paper

168 Scientia et Technica Año XXVIII, Vol. 29, No. 04, octubre-diciembre de 2024. Universidad Tecnológica de Pereira.

transformada wavelet [10], el uso de transporte hibrido [11], y

la reconstrucción de la reactividad por medio de detectores de

neutrones [12]. Un trabajo más reciente con mayor precisión

en el cálculo de la reactividad se presenta [13], allí se

discretiza la integral que contiene la dependencia de la

densidad de neutrones en la ecuación inversa de la cinética

puntual con infinitos números de Bernoulli.

En este trabajo se discretiza la dependencia de la densidad

de neutrones mediante la fórmula de Euler-Maclaurin y se

hace una aproximación únicamente con los tres primeros

números de Bernoulli resultado en una expresión donde la

corrección tiene características especiales como es el carácter

recursivo.

II. ASPECTOS TEÓRICOS

La seguridad en la generación de energía a partir de la fisión

nuclear es de vital importancia por lo que los modelos

computacionales que se utilicen para medir, modelar y/o

simular el comportamiento de un reactor nuclear deben

aproximarse en gran medida al comportamiento real. Para

esto, se parte de las ecuaciones de la cinética puntual que son

un conjunto de siete ecuaciones fuertemente acopladas que

describen la evolución temporal de la densidad de neutrones.

Están descritas por las siguientes ecuaciones [14]:

( ) ( )

i i

dP t t

P t C t

r -b

é ù

= +

ê ú

ë û

(1)

( )

( ); 1,2,...,6

i i

dC t

C t P t i

+ ( ) = =

(1)

donde, P es la densidad del neutrón, C

es la concentración del

i-ésimo grupo de neutrones atrasados, ρ es la reactividad, Λ es

el tiempo de generación de neutrones instantáneos, β

es la

fracción efectiva del i-ésimo grupo de neutrones atrasados, β

es la fracción total de neutrones atrasados, λ

es la constante de

decaimiento del i-ésimo grupo de neutrones atrasados.

Las ecuaciones (1) y (2) están sujetas a unas condiciones

iniciales, las cuales hacen que para un tiempo inicial t=0 la

reactividad sea nula, esto indica que el reactor está en estado

de criticidad. Al resolver (2) por cualquier método de solución

de ecuaciones diferenciales, y despejar la reactividad en (1),

considerando las condiciones iniciales mencionadas

anteriormente, se obtiene la expresión de la ecuación inversa

de la cinética puntual.

( )

( ) (0)

( )

( ) ( )

( )

t t

i i

dP t P

t e

P t dt P t

e P t dt

P t

r b b

l b

- -

= + -

(2)

La ecuación (3) tiene dependencia de la densidad de la

población de neutrones, lo cual hace que su implementación

computacionalmente no sea efectiva para cálculos en tiempo

real. Es la forma como se implementan algoritmos para

construir medidores de reactividad digital.

En la siguiente sección se describe un método que

discretiza la dependencia del histórico de la densidad de

neutrones lo cual aumenta la precisión en el cálculo de la

reactividad y facilita su implementación computacional.

III. MÉTODO PROPUESTO

Para discretizar la dependencia de la integral en la ecuación

(3), se emplea una aproximación de la integral a una suma

mediante la fórmula de Euler-Maclaurin [15],

 

(2 1) (2 1)

( ) [ ] [0] [ ]

[ ] [0]

(2 )!

x x

F y dy F s F F n

F n F

- -

= + +

é ù

- -

ë û

(3)

La ecuación (4) permite hacer un cambio en la dependencia

de funciones continuas F(y) a funciones discretas F[n].

Si x=3, la ecuación (4) se reescribe de la siguiente forma:

 

(1) (1)

(3) (3)

(5) (5)

( ) [ ] [0] [ ]

[ ] [0]

F y dy F s F F n

F n F

= + +

é ù

- -

ë û

é ù

+ -

ë û

é ù

- -

ë û

(4)

donde B

=1/6, B

=1/30 y B

=1/42 representan los tres

primeros números de Bernoulli. La reactividad con la

aproximación de tres números de Bernoulli puede ser

calculada, sustituyendo la ecuación (5) en (3) resulta:

1 2 3

[ ]

Trap B B B

r r r r r

= + ++

(6)

donde ρ

Trap

representa la reactividad con la corrección de la

regla trapezoidal [17], esto es:

 

(1)

1 1

[ ] [ ]

[ ] [ ] [ ] [0] [0] [ ]

[ ] 2

Trap i

i i i

i s

n P n e

P n P n

h n s P s h n P h P n

P n

r b b

= =

= + -

é ù

ê ú

- - - +

ê ú

ë û

å å

(7)

, ρ

y ρ

representan las correcciones con el primero, el

segundo y el tercer número de Bernoulli, dados por:

(1) (1)

[0] [ ] [0] [ ]

12 [ ]

[ ] [0] [ ] [0]

i i

h P n h P n

P n

h n P h n P

é ù

ê ú

- -

ë û

(8)

( )

(3) (2) (1)

(1) (2) (3)

(3) (2) (1)

(1) (2) (3)

[0] [ ] 3 [0] [ ]

3 [0] [ ] [0] [ ]

720 [ ]

[ ] [0] 3 [ ] [0]

3 [ ] [0] [ ] [0]

i i

h P n h P n

P n

h n P h n P

é ù

ê ú

+ +

ê ú

= -

ê ú

- -

ê ú

- -

ê ú

ë û

(9)

Scientia et Technica Año XXVIII, Vol. 29, No. 04, octubre-diciembre de 2024. Universidad Tecnológica de Pereira. 169

(2 1) ( )

2 1

[0] [ ]

30240 [ ]

2 1

[ ] [0]

p k k

h P n

P n

h n P

- -

é ù

æ ö

ê ú

ç ÷

ê ú

è ø

ë û

ê ú

é ù

æ ö

ê ú

ç ÷

ê ú

è ø

ë û

(10)

Es posible escribir una expresión para el cálculo de la

reactividad en una forma más compacta en términos de

expresiones pares e impares obteniendo la siguiente ecuación

para la reactividad:

6 3

1 1

(2 2 1) (2 )

(2 2 ) (2 1)

(2 2 1) (2 )

[ ] [ ] ( 1) *

[ ] (2 )!

2 1

[0] [ ]

2 1

[0] [ ]

2 1

[ ] [0]

2 1

Trap

i p

p k k

T B

n n

P n p

h P n

h n P

r r

= =

- -

= + -

é ù

æ ö

ê ú

ç ÷

ê ú

è ø

ê ú

æ ö

ê ú

ç ÷

ê ú

è ø

ë û

æ ö

ç ÷

è ø

å å

(2 2 ) (2 1)

[ ] [0]

p k k

h n P

- -

é ù

ê ú

é ù

ê ú

ë û

(11)

Se puede simplificar las derivadas de la respuesta del

sistema h

para un impulso unitario contenidas en la ecuación

(11) por ser una función exponencial de la siguiente forma:

( ')

(1) 2

( ')

(2) 3 2

( ')

(3) 4 3

( ') ( ')

t t

i i i i i

t t

i i i i i

t t

i i i i i

h t t e h t t

l b l

- -

- = = -

(12)

De otra parte, en la ecuación (11) se tienen términos pares e

impares de la densidad de neutrones P, para simplificar se

hace uso de la regla de derivadas para términos pares como

impares que fueron obtenidas integrando por partes la integral

de la ecuación (3), originando una serie geométrica que

converge [13], las cuales son válidas para diferentes formas de

la densidad de la población de neutrones:

(2)

(2 ) 2

( )

( ) ( ) ( )

( )

k k

P t

P t P t P t a

P t

é ù

= º

ê ú

ë û

(13)

(2)

(2 1) (1) (1) 2( 1)

( )

( ) ( ) ( )

( )

k k

P t

P t P t P t a

P t

- -

é ù

= º

ê ú

ë û

(14)

siendo a una constante que es la razón de convergencia de una

serie geométrica dada por:

(2)

( )

P t

a constante

P t

é ù

= =

ê ú

ë û

(15)

Con las consideraciones dadas por las ecuaciones (14-15) y

usando las ecuaciones (11-12) se puede escribir la reactividad

como:

6 3

1 1

2 2 1 2

(1)

2 2 2 1

2 2 1 2

(1)

[ ] [ ] ( 1) *

[ ] (2 )!

2 1

[ ]

[0]

2 1

[ ]

2 1

[0]

[ ]

2 1

[0]

Trap

i p

p k k

T B

n n

P n p

P n a

P n

P a

h n

r r

= =

- -

= + -

é ù

æ ö

ê ú

ç ÷

ê ú

è ø

ê ú

æ ö

ê ú

ç ÷

ê ú

è ø

ë û

æ ö

ç ÷

è ø

å å

2 2 2 1

p k k

- -

é ù

ê ú

é ù

ê ú

æ ö

ê ú

ç ÷

ê ú

è ø

ê ú

ë û

(16)

La ecuación (16) representa el método para calcular la

reactividad en forma discreta con la dependencia de solo hasta

la segunda derivada, siendo el orden de las derivadas menor al

presentado en las ecuaciones (8-10).

Una representación más simple y compacta para escribir la

ecuación (16) es:

 

[ ]

Trap C

n n

r r r

= +

(17)

La ecuación (17) tiene dos ventajas, la primera es el carácter

semi recursiva, lo que no se obtenía con la ecuación inversa de

la cinética puntual dada por la ecuación (3) ya que es una

integral de convolución y depende de todos los puntos de la

densidad de neutrones. La segunda es que los términos

𝜌

𝐶

𝑃

que originalmente fueron definidos en las ecuaciones (8-10),

ya no dependen de las derivadas de orden superiores, sino que

su nueva dependencia solo es función de hasta la segunda

derivada de la densidad de la población de neutrones.

IV. RESULTADOS

En esta sección se exponen los resultados de las

simulaciones numéricas realizadas a partir de la ecuación (17).

Las simulaciones se presentan para un reactor nuclear que

utiliza combustible

235

𝑈

(los parámetros de este tipo de

reactor están dados en la Tabla I).

TABLA I. PARÁMETROS DEL COMBUSTIBLE

235

UTILIZADOS EN LOS

EXPERIMENTOS NUMÉRICOS

Parámetro

Valor

Parámetro

Valor [s

-1

]

𝛽

0.000266

𝜆

0.0127

𝛽

0.001491

𝜆

0.0317

𝛽

0.001316

𝜆

0.115

𝛽

0.002849

𝜆

0.311

𝛽

0.000896

𝜆

1.4

𝛽

0.000182

𝜆

3.87

El tiempo de generación de neutrones es Λ = 2 ×10

-5

s. Los

resultados se comparan con la solución exacta o método de

170 Scientia et Technica Año XXVIII, Vol. 29, No. 04, octubre-diciembre de 2024. Universidad Tecnológica de Pereira.

referencia usando la ecuación (3), también se realizan

comparaciones con la aproximación del primer y segundo

número de Bernoulli. y otros métodos encontrados en la

literatura como el método de Hamming [9] y el método de

Lagrange 5 puntos [7]. Las simulaciones se realizan para

diferentes densidades de neutrones y diferentes pasos de

tiempo.

El primer experimento numérico se realiza para una

densidad de neutrones de la forma P(t)=exp(ωt) para

diferentes pasos de tiempo, el valor de la frecuencia

considerada ω y el tiempo de simulación son: ω=0.0243 s

-1

t=1000 s respectivamente. Los resultados en las máximas

diferencias se presentan en la Tabla II. Es posible observar que

para los valores ω, el método propuesto usando (6), esto es el

cálculo de la reactividad incluyendo todas las derivadas, la

aproximación del tercer número de Bernoulli representado por

𝜌

𝐵

, tiene mayor precisión que los dos primeros números de

Bernoulli denotados por

𝜌

𝐵

𝜌

𝐵

. Iguales resultados son

encontrados cuando se compara con el método reduciendo el

orden de las derivadas

𝜌

𝐶

usando (17), estos resultados son

equivalentes pues la forma de la densidad de neutrones cumple

con las condiciones dadas por las ecuaciones (13-15).

TABLA II. MÁXIMA DIFERENCIA EN LA REACTIVIDAD EN PCM

En las Tablas III-IV se presentan los resultados numéricos

obtenidos para frecuencias y tiempos de simulación:

ω=0.01046 s

-1

y t = 800 s; ω=0.02817 s

-1

y t=600 s

respectivamente. Se observa que para los valores diferentes ω,

el método propuesto usando la aproximación del tercer

número de Bernoulli representado por

𝜌

𝐵

, tiene mayor

precisión que los métodos de Lagrange de 5 puntos y el

método de Hamming, excepto para un paso de tiempo de

cálculo muy alto de T = 1s. Iguales resultados son encontrados

cuando se compara con el método reduciendo el orden de las

derivadas

𝜌

𝐶

ya que cumple con las condiciones dadas por las

ecuaciones (13-15).

TABLA III. MÁXIMA DIFERENCIA EN LA REACTIVIDAD EN PCM

T(s)

Lagrange 5 puntos

Hamming

𝜌

𝐵

𝜌

𝐶

0.01

2.37×10

-5

1.062×10

-10

3.77×10

-11

3.77×10

-11

0.05

3.98×10

-4

2.43×10

-7

4.46×10

-11

4.46×10

-11

0.1

1.14×10

-3

5.29×10

-6

7.71 ×10

-9

7.71×10

-9

0.5

1.04×10

-2

7.86×10

228

2.80×10

-3

2.80×10

-3

1.83×10

-2

Infinito

5.60×10

-1

5.60×10

-1

TABLA IV. MÁXIMA DIFERENCIA EN LA REACTIVIDAD EN PCM

T(s)

Lagrange 5 puntos

Hamming

𝜌

𝐵

𝜌

𝐶

0.01

6.38×10

-5

2.88×10

-10

2.96×10

-11

2.96×10

-11

0.05

1.06×10

-3

6.56×10

-07

4.10×10

-11

4.01×10

-11

0.1

3.05×10

-3

1.42×10

-05

7.95×10

-9

7.95×10

-9

0.5

2.65×10

-2

4.44×10

182

2.80×10

-3

2.80×10

-3

4.36×10

-2

Infinito

5.77×10

-1

5.77×10

-1

Es bueno resaltar que, para cualquier paso de cálculo, el

método propuesto ya sea reduciendo o no el orden de las

derivadas, él método propuesto tiene mejores aproximaciones

que el método de Hamming con un paso de cálculo de T=1 s,

el método de Hamming diverge, mientras los métodos

𝜌

𝐵

𝜌

𝐶

convergen a valores cercanos de 10

-1

pcm. Esto es posible

debido a la alta precisión que proporciona la aproximación del

tercer número de Bernoulli en el cálculo de la reactividad.

La Figura 1 muestra la curva de reactividad en pcm para

una densidad de neutrones de forma cúbica P(t)=a+bt

, a=1 y

b= (0.0127)5/9 y un paso de cálculo T=1 s para un tiempo

total de simulación t=10000 s. Se observa que el tercer

número de Bernoulli denotado por

𝜌

𝐵

tiene mejores

aproximaciones cuando se comparan con el método de

referencia dado por la solución exacta que se obtiene al

resolver analíticamente (3).

0 2000 4000 6000 8000 10000

-2

Reactividad [pcm]

Tiempo [s]

Método de referencia

Fig. 1. Curva de reactividad en pcm, P(t)=a+bt

, a=1, b= (0.0127)

En los próximos experimentos numéricos se considera un

tiempo total de simulación de t=10000 s, en las Tablas V-VI

se exponen los resultados obtenidos cuando se usan las

derivadas normales y cuando se reduce el orden

respectivamente para una densidad de neutrones de la forma

P(t)=a+sin(bt) con distintos valores de a y fijando b= π/10.

Esta forma de la densidad de neutrones no cumple

estrictamente con las condiciones dadas por (13-15). Es

posible notar que hay cambios en los resultados usando el

tercer número de Bernoulli,

𝜌

𝐵

𝜌

𝐶

cuando se reduce el

T(s)

𝜌

𝐵

𝜌

𝐵

𝜌

𝐵

𝜌

𝐶

0.01

6.17×10

-8

4.50×10

-11

4.30×10

-11

4.30×10

-11

0.05

3.85×10

-5

3.20×10

-8

4.69×10

-11

4.69×10

-11

0.1

6.14×10

-4

2.04×10

-6

7.60×10

-9

7.60×10

-9

0.2

9.73×10

-3

1.29×10

-4

1.92×10

-6

1.92×10

-6

0.5

3.56×10

-1

2.93×10

-2

2.70×10

-3

2.70×10

-3

4.67×10

1.50×10

5.53×10

-1

5.53×10

-1

Scientia et Technica Año XXVIII, Vol. 29, No. 04, octubre-diciembre de 2024. Universidad Tecnológica de Pereira. 171

orden de la derivada hasta 5. Las grandes diferencias aparecen

cuando se considera el valor a=1 ya que los resultados

obtenidos pasan de estar en el orden de 10

-6

pcm cuando no se

reduce el orden de las derivadas y cuando se reducen las

máximas diferencias alcanzan valores de 10

-1

pcm.

TABLA V. MÁXIMA DIFERENCIA EN LA REACTIVIDAD EN PCM

𝑇

(s)

𝜌

𝐵

𝜌

𝐵

𝜌

𝐵

0.1

5.09×10

-2

3.44×10

-4

2.51×10

-6

0.1

100

6.15×10

-4

2.05×10

-6

7.63×10

-9

0.1

150

6.14×10

-4

2.05×10

-6

7.62×10

-9

TABLA VI. MÁXIMA DIFERENCIA EN LA REACTIVIDAD EN PCM

𝑇

(s)

𝜌

𝐶

𝜌

𝐶

𝜌

𝐶

0.1

5.09×10

-2

4.04×10

-1

3.83×10

-1

0.1

100

6.15×10

-4

2.06×10

-6

2.06×10

-8

0.1

150

6.14×10

-4

2.05×10

-6

1.66×10

-8

En el último experimento numérico se puede observar el

comportamiento del método propuesto cuando se usan la

reducción del orden de las derivadas para una forma de la

densidad de la población de neutrones P(t)=a+sin(bt) con

a=150, b= π/10 y un paso de tiempo T=1 s. Las máximas

diferencias en reactividad son aproximadamente 4.67 pcm,

1.50 pcm y 0.55 pcm usando el primero, segundo y tercer

número de Bernoulli respectivamente. Los resultados indican

nuevamente como el caso

𝜌

𝐶

presenta buenas

aproximaciones.

En la Figura 2 se muestra la curva de reactividad para una

densidad de neutrones P(t)=a +b sinh(ka t) con a=100, b= 1 y

ka= 1.27×10

-3

-1

con un tamaño de paso de tiempo T=1 s

hasta un tiempo de simulación de t=10000 s. Se evidencia que

el método propuesto tiene una mejor aproximación en

comparación con el primer y el segundo número de Bernoulli.

0 2000 4000 6000 8000 10000

-2

Reactividad [pcm]

Tiempo [s]

Método de referencia

Fig. 2. Curva de reactividad en pcm para P(t)=a +bsinh(ka t), t=10000 s.

V. CONCLUSIONES

Se presento un método para calcular numéricamente la

reactividad considerando los parámetros físicos para un

reactor nuclear que emplea combustible

235

𝑈

usando la

ecuación inversa de la cinética puntual, la cual depende de la

densidad de población de neutrones; esto se logró empleando

la fórmula de Euler-Maclaurin que permite hacer una

aproximación de una integral entre el tiempo continuo y el

tiempo discreto. El método propuesto usó la aproximación del

tercer número de Bernoulli y el criterio de reducir el orden de

las derivadas, permitiendo calcular la reactividad en una forma

semi recursiva, obteniendo buenas aproximaciones para

diferentes formas de la densidad de la población de neutrones

y diferentes pasos de cálculo. La limitación y las dificultades

se pueden presentar por las aproximaciones dadas en las

ecuaciones (13-15). Sin embargo, existe buena precisión y se

evidencia en los resultados numéricos realizados, la

aproximación del tercer número de Bernoulli reduce el costo

computacional y puede ser un método alternativo para

implementar en medidores de reactividad digital.

REFERENCIAS

[1] M. Stacey, Nuclear Reactor Physics 3e. Weinheim,

Germany: Wiley-VCH Verlag GmbH & Co. KGaA,

2018. doi: 10.1002/9783527812318

[2] Y. Shimazu, Y. Nakano, Y. Tahara, and T. Okayama,

“Development of a Compact Digital Reactivity Meter and

a Reactor Physics Data Processor.,” Nucl Technol, vol.

77, no. 3, pp. 247–254, 1987, doi: 10.13182/NT87-

A33964

[3] J. E. Hoogenboom and A. R. van der Sluijs, “Neutron

source strength determination for on-line reactivity

measurements,” Ann Nucl Energy, vol. 15, no. 12, pp.

553–559, 1988, doi: 10.1016/0306-4549(88)90059-X

[4] S. E. Binney and A. J. M. Bakir, “Design and

development of a personal-computer-based reactivity

meter for a research reactor,” Nucl Technol, vol. 85, no.

1, pp. 12–21, 1989, doi: 10.13182/NT89-A34223

[5] S. A. Ansari, “Development of on-line reactivity meter

for nuclear reactors,” IEEE Trans Nucl Sci, vol. 38, no. 4,

pp. 946–952, Aug. 1991, doi: 10.1109/23.83857

[6] A. Kitano, M. Itagaki, and M. Narita, “Memorial-Index-

Based Inverse Kinetics Method for Continuous

Measurement of Reactivity and Source Strength,” J Nucl

Sci Technol, vol. 37, no. 1, pp. 53–59, Jan. 2000, doi:

10.1080/18811248.2000.9714866

[7] H. Malmir and N. Vosoughi, “On-line reactivity

calculation using Lagrange method,” Ann Nucl Energy,

vol. 62, pp. 463–467, 2013, doi:

10.1016/j.anucene.2013.07.006

[8] D. Suescún-Díaz, J. A. Rodríguez-Sarasty, and J. H.

Figueroa-Jiménez, “Reactivity calculation using the

Euler–Maclaurin formula,” Ann Nucl Energy, vol. 53, pp.

104–108, Mar. 2013, doi: 10.1016/j.anucene.2012.09.026

[9] D. Suescún-Díaz, M. C. Ibarguen-Gonzalez, and J. H.

Figueroa-Jiménez, “Hamming generalized corrector for

172 Scientia et Technica Año XXVIII, Vol. 29, No. 04, octubre-diciembre de 2024. Universidad Tecnológica de Pereira.

reactivity calculation,” Kerntechnik, vol. 79, no. 3, pp.

219–225, Jun. 2014, doi: 10.3139/124.110423

[10] M. Mohideen Abdul Razak and N. Rathinasamy, “Haar

wavelet for solving the inverse point kinetics equations and

estimation of feedback reactivity coefficient under

background noise,” Nuclear Engineering and Design, vol.

335, pp. 202–209, Aug. 2018, doi:

10.1016/j.nucengdes.2018.04.022

[11] P. Picca and R. Furfaro, “Reactivity determination using

the hybrid transport point kinetics and the area method,”

Ann Nucl Energy, vol. 114, pp. 191–197, 2018, doi:

10.1016/j.anucene.2017.12.019

[12] N. Chentre, P. Saracco, S. Dulla, and P. Ravetto, “On the

prompt time eigenvalue estimation for subcritical

multiplying systems,” Ann Nucl Energy, vol. 132, pp. 172–

180, 2019, doi: 10.1016/j.anucene.2019.04.030

[13] D. Suescún-Díaz, G. Ule-Duque, and F. H. Escobar,

“Novel approach to solving the inverse equation of point

kinetics by the Bernoulli number generalisation method,” J

Nucl Sci Technol, vol. 57, no. 8, pp. 989–999, 2020, doi:

10.1080/00223131.2020.1742813

[14] J. J. Duderstadt and L. J. Hamilton, Nuclear Reactor

Analysis, Second ed. New York: John Wiley & Sons Inc,

1976.

[15] Y. K. Kwok, Applied Complex Variables for Scientists and

Engineers. Cambridge University Press, 2010. doi:

10.1017/CBO9780511844690

[16] S. Haykin and B. Van Veen, Signals and Systems. New

York: Wiley, 2005.

[17] D. Suescún Díaz, A. Senra Martinez, and F. Carvalho Da

Silva, “Calculation of reactivity using a finite impulse

response filter,” Ann Nucl Energy, vol. 35, no. 3, pp. 472–

477, 2008, doi: 10.1016/j.anucene.2007.07.002

Daniel Suescún Díaz. Recibió el título de Licenciando en

Matemáticas de la Universidad Industrial de Santander en

1993, el título de Físico de la Universidad Industrial de

Santander en 1999 y MSc en Física de la Universidad

Industrial de Santander en el año 2000, el título de Doctor en

Física en Reactores Nucleares de la Universidad Federal de

Rio de Janeiro en el año 2007. Entre sus intereses

investigativos se encuentra la Física nuclear y la Física

computacional con métodos numéricos estocásticos.

https://orcid.org/0000-0003-2422-0684

Geraldyne Ule Duque. Recibió el título de Físico de la

Universidad Surcolombiana en 2019, el título M Sc. en Física

de la Universidad de Valencia en el año 2021, actualmente

realiza estudios PhD en Física en la Universidad de Bangor,

Reino Unido. Entre los intereses investigativos se encuentra la

Física computacional con métodos numéricos, método de

Monte Carlo, procesos estocásticos, Física de reactores

nucleares.

https://orcid.org/0000-0001-9534-5260

Diego Peña Lara. Recibió el título de Físico de la Universidad

del Valle 1986, el título M Sc. en Física de la Universidad de.

Valle en el año 1990 y el título de Doctor en Física de la

Universidad Federal de Minas Gerais en el año 1999. Entre los

intereses investigativos se encuentra la Física Estadística,

Física computacional con métodos numéricos, dinámica

molecular, método de Monte Carlo, procesos estocásticos,

Transiciones de fases en sistemas magnético e iónicos

https://orcid.org/0000-0001-6199-1547