from mpi4py import MPI
import numpy as np
def trapezoidal_rule(func, a, b, n):
h = (b - a) / n
integral = (func(a) + func(b)) / 2.0
for i in range(1, n):
x_i = a + i * h
integral += func(x_i)
integral *= h
return integral
def midpoint_rule(func, a, b, n):
h = (b - a) / n
integral = 0.0
for i in range(n):
x_mid = a + (i + 0.5) * h
integral += func(x_mid)
integral *= h
return integral
if __name__ == "__main__":
comm = MPI.COMM_WORLD
rank = comm.Get_rank()
size = comm.Get_size()
if rank == 0:
# User input
function_str = input("Enter the function (in terms of x): ")
a = float(input("Enter the lower bound (a): "))
b = float(input("Enter the upper bound (b): "))
n_values = [int(input(f"Enter the number of intervals (n) for experiment {i+1}: ")) for i in range(3)]
else:
function_str = None
a = None
b = None
n_values = [None] * 3
# Broadcast user input to all processes
function_str = comm.bcast(function_str, root=0)
a = comm.bcast(a, root=0)
b = comm.bcast(b, root=0)
n_values = comm.bcast(n_values, root=0)
# Parse the function string and create a function object
def func_wrapper(x):
return eval(function_str)
for n in n_values:
if rank == 0:
print(f"Experiment with n = {n}:")
local_n = n // size
local_a = a + rank * (b - a) / size
local_b = local_a + (b - a) / size
local_integral_trapezoidal = trapezoidal_rule(func_wrapper, local_a, local_b, local_n)
local_integral_midpoint = midpoint_rule(func_wrapper, local_a, local_b, local_n)
integral_trapezoidal = comm.reduce(local_integral_trapezoidal, op=MPI.SUM, root=0)
integral_midpoint = comm.reduce(local_integral_midpoint, op=MPI.SUM, root=0)
if rank == 0:
print("Trapezoidal Rule:")
print(f" Computed Integral: {integral_trapezoidal}")
print(f" Time taken: {end_time - start_time} seconds")
print("Midpoint Rule:")
print(f" Computed Integral: {integral_midpoint}")
print(f" Time taken: {end_time - start_time} seconds")
MPI.Finalize()
ZnJvbSBtcGk0cHkgaW1wb3J0IE1QSQppbXBvcnQgbnVtcHkgYXMgbnAKaW1wb3J0IHRpbWUKCmRlZiB0cmFwZXpvaWRhbF9ydWxlKGZ1bmMsIGEsIGIsIG4pOgogICAgaCA9IChiIC0gYSkgLyBuCiAgICBpbnRlZ3JhbCA9IChmdW5jKGEpICsgZnVuYyhiKSkgLyAyLjAKICAgIGZvciBpIGluIHJhbmdlKDEsIG4pOgogICAgICAgIHhfaSA9IGEgKyBpICogaAogICAgICAgIGludGVncmFsICs9IGZ1bmMoeF9pKQogICAgaW50ZWdyYWwgKj0gaAogICAgcmV0dXJuIGludGVncmFsCgpkZWYgbWlkcG9pbnRfcnVsZShmdW5jLCBhLCBiLCBuKToKICAgIGggPSAoYiAtIGEpIC8gbgogICAgaW50ZWdyYWwgPSAwLjAKICAgIGZvciBpIGluIHJhbmdlKG4pOgogICAgICAgIHhfbWlkID0gYSArIChpICsgMC41KSAqIGgKICAgICAgICBpbnRlZ3JhbCArPSBmdW5jKHhfbWlkKQogICAgaW50ZWdyYWwgKj0gaAogICAgcmV0dXJuIGludGVncmFsCgppZiBfX25hbWVfXyA9PSAiX19tYWluX18iOgogICAgY29tbSA9IE1QSS5DT01NX1dPUkxECiAgICByYW5rID0gY29tbS5HZXRfcmFuaygpCiAgICBzaXplID0gY29tbS5HZXRfc2l6ZSgpCgogICAgaWYgcmFuayA9PSAwOgogICAgICAgICMgVXNlciBpbnB1dAogICAgICAgIGZ1bmN0aW9uX3N0ciA9IGlucHV0KCJFbnRlciB0aGUgZnVuY3Rpb24gKGluIHRlcm1zIG9mIHgpOiAiKQogICAgICAgIGEgPSBmbG9hdChpbnB1dCgiRW50ZXIgdGhlIGxvd2VyIGJvdW5kIChhKTogIikpCiAgICAgICAgYiA9IGZsb2F0KGlucHV0KCJFbnRlciB0aGUgdXBwZXIgYm91bmQgKGIpOiAiKSkKICAgICAgICBuX3ZhbHVlcyA9IFtpbnQoaW5wdXQoZiJFbnRlciB0aGUgbnVtYmVyIG9mIGludGVydmFscyAobikgZm9yIGV4cGVyaW1lbnQge2krMX06ICIpKSBmb3IgaSBpbiByYW5nZSgzKV0KICAgIGVsc2U6CiAgICAgICAgZnVuY3Rpb25fc3RyID0gTm9uZQogICAgICAgIGEgPSBOb25lCiAgICAgICAgYiA9IE5vbmUKICAgICAgICBuX3ZhbHVlcyA9IFtOb25lXSAqIDMKCiAgICAjIEJyb2FkY2FzdCB1c2VyIGlucHV0IHRvIGFsbCBwcm9jZXNzZXMKICAgIGZ1bmN0aW9uX3N0ciA9IGNvbW0uYmNhc3QoZnVuY3Rpb25fc3RyLCByb290PTApCiAgICBhID0gY29tbS5iY2FzdChhLCByb290PTApCiAgICBiID0gY29tbS5iY2FzdChiLCByb290PTApCiAgICBuX3ZhbHVlcyA9IGNvbW0uYmNhc3Qobl92YWx1ZXMsIHJvb3Q9MCkKCiAgICAjIFBhcnNlIHRoZSBmdW5jdGlvbiBzdHJpbmcgYW5kIGNyZWF0ZSBhIGZ1bmN0aW9uIG9iamVjdAogICAgZGVmIGZ1bmNfd3JhcHBlcih4KToKICAgICAgICByZXR1cm4gZXZhbChmdW5jdGlvbl9zdHIpCgogICAgc3RhcnRfdGltZSA9IHRpbWUudGltZSgpCgogICAgZm9yIG4gaW4gbl92YWx1ZXM6CiAgICAgICAgaWYgcmFuayA9PSAwOgogICAgICAgICAgICBwcmludChmIkV4cGVyaW1lbnQgd2l0aCBuID0ge259OiIpCiAgICAgICAgbG9jYWxfbiA9IG4gLy8gc2l6ZQogICAgICAgIGxvY2FsX2EgPSBhICsgcmFuayAqIChiIC0gYSkgLyBzaXplCiAgICAgICAgbG9jYWxfYiA9IGxvY2FsX2EgKyAoYiAtIGEpIC8gc2l6ZQoKICAgICAgICBsb2NhbF9pbnRlZ3JhbF90cmFwZXpvaWRhbCA9IHRyYXBlem9pZGFsX3J1bGUoZnVuY193cmFwcGVyLCBsb2NhbF9hLCBsb2NhbF9iLCBsb2NhbF9uKQogICAgICAgIGxvY2FsX2ludGVncmFsX21pZHBvaW50ID0gbWlkcG9pbnRfcnVsZShmdW5jX3dyYXBwZXIsIGxvY2FsX2EsIGxvY2FsX2IsIGxvY2FsX24pCgogICAgICAgIGludGVncmFsX3RyYXBlem9pZGFsID0gY29tbS5yZWR1Y2UobG9jYWxfaW50ZWdyYWxfdHJhcGV6b2lkYWwsIG9wPU1QSS5TVU0sIHJvb3Q9MCkKICAgICAgICBpbnRlZ3JhbF9taWRwb2ludCA9IGNvbW0ucmVkdWNlKGxvY2FsX2ludGVncmFsX21pZHBvaW50LCBvcD1NUEkuU1VNLCByb290PTApCgogICAgICAgIGlmIHJhbmsgPT0gMDoKICAgICAgICAgICAgZW5kX3RpbWUgPSB0aW1lLnRpbWUoKQogICAgICAgICAgICBwcmludCgiVHJhcGV6b2lkYWwgUnVsZToiKQogICAgICAgICAgICBwcmludChmIiAgQ29tcHV0ZWQgSW50ZWdyYWw6IHtpbnRlZ3JhbF90cmFwZXpvaWRhbH0iKQogICAgICAgICAgICBwcmludChmIiAgVGltZSB0YWtlbjoge2VuZF90aW1lIC0gc3RhcnRfdGltZX0gc2Vjb25kcyIpCiAgICAgICAgICAgIHByaW50KCJNaWRwb2ludCBSdWxlOiIpCiAgICAgICAgICAgIHByaW50KGYiICBDb21wdXRlZCBJbnRlZ3JhbDoge2ludGVncmFsX21pZHBvaW50fSIpCiAgICAgICAgICAgIHByaW50KGYiICBUaW1lIHRha2VuOiB7ZW5kX3RpbWUgLSBzdGFydF90aW1lfSBzZWNvbmRzIikKCiAgICBNUEkuRmluYWxpemUoKQ==