how to find an array made by results of integrals-CodePudding

In my code res1 is the result of an integral. My goal is to repeat this integral for 1000 times (with for loop maybe but I can not figure out how) with one different variable which is :

Energy = np.arange(2.1,102.1,0.1)

so in the code below, E means energy but it is the first variable of that array. I want to have an array that is made of the results of the integrals ( by changin E like 2.1 to 2.2 and then to 2.3)

import numpy as np
from math import pi,sqrt
vtheta=np.linspace(0.0,pi,1000)
def f(theta):
    E=2.1
    val=(2.0)/np.sqrt(E-1 np.cos(theta))
    result=(val.sum()-(sqrt(1./2.1) sqrt(1./298.)))*((np.pi-0)/999.)
    return result
res1=f(vtheta)

when I try to put the for loop in the function I change the inside of result which I do not want to do.

CodePudding user response：

Something like this?

import numpy as np
import matplotlib.pyplot as plt

def f(theta, E=2.1):
    val = 2. / np.sqrt(E - 1   np.cos(theta))
    return (val.sum() - np.sqrt(1. / 2.1) - np.sqrt(1. /298.)) * (np.pi / 999.)

vtheta = np.linspace(0., np.pi, 1000)
energies = np.arange(2.1, 102.1, .1)
results = []
for e in energies:
    results.append(f(vtheta, E=e))
    
plt.plot(energies, results)
plt.xlabel('Energy')
plt.show()

Gives you the following plot:

CodePudding user response：

Try this, with two different implementations, u can use either f1() or f2()

import numpy as np
from math import pi, sqrt
import time

vtheta = np.linspace(0.0, pi, 1000)
E = np.arange(2.1, 102.1, 0.1)


def f1(theta,E):

    # using list comprehension
    val = [(2.0 / np.sqrt(i - 1   np.cos(theta))) for i in E]
    result = [((i.sum() - (sqrt(1. / 2.1)   sqrt(1. / 298.))) * ((np.pi - 0) / 999.)) for i in val]

    return result

def f2(theta,E):

    # using looping
    result = []

    for i in E:
        val = 2.0 / np.sqrt(i - 1   np.cos(theta))
        result.append((val.sum() - (sqrt(1. / 2.1)   sqrt(1. / 298.))) * ((np.pi - 0) / 999.))

    return result

if __name__ == "__main__":

    t0 = time.time()
    res1 = f1(vtheta,E)
    tf = time.time()
    print(f"Process time for list comprehension: {tf-t0}")

    t0 = time.time()
    res2 = f2(vtheta,E)
    tf = time.time()
    print(f"Process time for looping: {tf-t0}")

    if res2==res1:
        print("Both result matched")