La función gamma
Factorial de un número
Recordaremos en primer lugar, el concepto de factorial de un número entero positivo n.
n!=n·(n-1)·(n-2)...3·2·1
MATLAB dispone de una función denominada factorial, lo calculamos también con la función prod.
>> 1*2*3*4*5 ans = 120 >> factorial(5) ans = 120 >> prod(1:5) ans = 120
La función gamma
Se define del siguiente modo
Donde z puede ser un número real o complejo.
Esta función es más complicada en la parte negativa de z, tiene asíntotas verticales para z=0, -1, -2, -3...
syms x; ezplot(gamma(x),[-4,5]), grid on
Representamos gráficamente la función gamma para z>0, llamando a la función del mismo nombre disponible en MATLAB.
>> fplot(@gamma,[0.1,4]) >> grid on >> xlabel('x') >> ylabel('\gamma(x)') >> title('Función gamma')
Calculamos el mínimo de la función aproximadamente.
>> x=1:0.01:2; >> y=gamma(x); >> yp=diff(y); >> idx=find(abs(yp)<0.0001); >> x(idx) ans = 1.4500 1.4600 >> y(idx) ans = 0.8857 0.8856
El mínimo se encuentra entre 1 y 2 en la posición x=1.46 y vale 0.8856.
Relaciones
Relacionamos Γ(z+1) con Γ(z)
Integrando por partes
Mediante esta relación calculamos el valor de la función gamma, para z=1,2,3,4,...
>> z=1:6; >> gamma(z) ans = 1 1 2 6 24 120
Calculamos la función Γ(z) para valores fraccionarios de z. Por ejemplo, para z=1/2
>> syms t; >> y=exp(-t)/sqrt(t); >> int(y,0,inf) ans =pi^(1/2)
Para calcular el valor de la integral hacemos el cambio t=y2.
Se puede consultar en la tablas el valor de esta integral, o resolverla del siguiente modo
Cambiamos de coordenadas rectangulares a polares. En coordendas polares el elemento diferencial de área dx·dy=(r·dθ)dr
Finalmente, utilizando la relación anterior Γ(z+1)=z·Γ(z)
>> n=6; >> num=2*(1:n)-1; >> prod(num)*sqrt(pi)/2^n ans = 287.8853 >> gamma(n+1/2) ans = 287.8853
Fórmula de Euler
>> z=1.3; >> gamma(z)*gamma(1-z) ans = -3.8832 >> pi/sin(pi*z) ans = -3.8832
La función beta
La función beta está relacionada con la función gamma
Periodo de las oscilaciones
Consideremos una partícula de masa m oscilando a lo largo del eje X bajo la acción de una fuerza conservativa cuya energía potencial tiene la forma
b es un número real positivo y n suele ser un entero, pero puede ser un número real positivo
- El caso n=1, se ha estudiado en la página titulada Oscilaciones de una partícula bajo la acción de una fuerza de módulo constante.
- El caso n=2, corresponde a un oscilador armonico.
- El caso n=4, se ha estudiado en la página titulada Oscilador no lineal
Representamos la energía potencial Ep(x) entre -1.2<x<1.2, tomando b=1, para cada uno de los siguientes casos: n=1,2,3,4,5. Se trata de una función simétrica, por lo que solamente es necesario calcular los valores entre 0 y 1.2
x=linspace(0,1.2,200); hold on for n=1:5 y=x.^n; yy=[fliplr(y),y]; xx=[-fliplr(x),x]; plot(xx,yy,'displayName',num2str(n)) end hold off legend('-DynamicLegend','location','best') grid on ylim([0,1.2]) xlabel('x') ylabel('E_p(x)') title('Energía potencial')
La energía constante E de la partícula es
Para calcular el periodo P de la oscilación despejamos dt e integramos
El periodo es cuatro veces el tiempo que tarda en desplazarse desde el origen hasta A, la amplitud de la oscilación, la posición en la que la velocidad dx/dt de la partícula es nula. E=bAn
En el intervalo 0≤x≤A, x es positivo por lo que |x|=x
Para resolver la integral hacemos el cambio de variable
Por otra parte,
cuando x=A, sin2θ=1, θ=π/2
En términos de la nueva variable θ, el periodo se expresa
Buscamos esta integral en una tabla de integrales, pág. 397, 3.621, n° 1
con μ=2/n. Donde B es la función beta, teniendo en cuenta su relación con la función gamma
Obtenemos la fórmula del periodo en términos de la energía total E de la partícula
En términos de la amplitud A
Representamos , en función de la amplitud A
hold on for n=1:4 f=@(x) 2^(2/n)*sqrt(2)*x.^(1-n/2)*gamma(1/n)^2/(n*gamma(2/n)); fplot(f,[0,5],'displayName',num2str(n)) end hold off grid on legend('-DynamicLegend','location','best') ylim([0,25]) xlabel('A') ylabel('P') title('Periodo en función de la amplitud')
Para n=2, el periodo es independiente de la amplitud
Ejemplos
Para n=1
Para n=2
Para n=4
En la página titulada Oscilaciones de una partícula bajo la acción de una fuerza de módulo constante, b=f·m
Para el sistema formado por una masa m unida a un muelle elástico de constante k, Ep(x)=kx2/2, por lo que b=k/2
La energía potencial correspondiente a la fuerza conservativa proporcional a x3 es
El periodo es
En la página titulada Oscilador no lineal, obtuvimos la fórmula del periodo
Comprobamos que dan el mismo resultado
>> 4*ellipke(1/2) ans = 7.4163 >> gamma(1/4)^2/sqrt(pi) ans = 7.4163
Ejemplos en el Curso de Física
Vaciado de un depósito abierto.
Referencias
Special functions and their applications. N.N. Lebedev. Prentice Hall Inc. 1965
Anyi Amezquita, María Delgado, Diego Rasero. Sistemas oscilantes intrínsecamente no lineales. Revista Brasileira de Ensino de Física, vol. 43, e20210110 (2021)
Pirooz Mohazzabi. Theory and examples of intrinsically nonlinear oscillators. Am. J. Phys. 72 (4), April 2004, pp. 492-498
I. S. Gradshteyn, I. M. Ryzhik. Table of Integrals, Series, and Products. Seventh Edition. Elsevier (2007). Pág. 397, 3.621, n° 1