Using summation to calculate for partition sums

Currently in the process of calculating for this formula in MATLAB:

A bit confused on how to set it up though. I started with symsum but had some issues understanding how the variables syms variable is declared.

For the values within the formula, it is a formula for calculating the total internal partition sum. The total internal partition sum is given by summing the factor e^−hcFs/kT over all states s of an isotopologue. Here hcFs is the total energy (electronic, vibrational, rotational, and any other quantized motion) of state s. di and ds are degeneracy factors that are state independent and state dependent, respectively. T is temperature. Here is the code I have so far.

clear all

wn = 3057.169717; % wavenumber of CH4 chosen
h = 6.62607015*10^-27; % Planck's constant
c = 2.99792458*10^10; % speed of light
k = 1.380649*10^-16; % Boltzmann constant
c2 = 1.4387769; % second radiation constant
En = 1780.1223; % lower state energy
d1 = 117; % upper state degeneracy
d2 = 111; % lower state degeneracy
syms d1 d2

tempF = 70;
tempK = ((5/9)*(tempF + 459.67)); % temp in Kelvin

Q(tempK) = (d1).*symsum((d2).*exp(-(En)/(k)*(tempK)), d1, d2)

I want to get an output number for Q(tempK), not a row of 0's.

[Torsten]

I wonder what state vectors you want to sum over. I see only scalars, no vectors in your code.

Answer 1

[Vinayak Choyyan]

Hi zack

As per my understanding, you are trying to sum over all the states of an isotopologue. ds is the degeneracy factor of the state dependent, and is probably an vector containing values for each state.

You tried to use symsum(), but this function increments the variable it is summing by 1, just like in regular mathematics summation operation. Also,

symsum(f,k,a,b)

returns the sum of the series f with respect to the summation index k from the lower bound a to the upper bound b. You have given d1, d2 but d1 is greater than d2. You can read more about symsum here.

Here is an example that might help you solve the equation you are trying to compute.

d1=sym(5);
d2=sym([1 2 3 4]);
En = 1780.1223;
k = 1.380649*10^-16;
tempF = 70;
tempK = ((5/9)*(tempF + 459.67));
Q = (d1)*sum(d2*exp(-(En)/(k*tempK)))

[zen]

Thank you very much, this has helped me get on the right track.