1a

from pylab import *

x=[1,3,5,7,9]

y=[2,4,6,8,10]

if (len(x)!=len(y)):

print("Values don't match")

else:

scatter(x,y)

plot(x,y)

xlabel("X-axis")

ylabel('Y-Axis')

title('Graph connecting the given plots')

axvline(x=0,color='r')

axhline(y=0,color='g')

grid()

show()

1 b.

from sympy import *

x,y= symbols('x,y')

f= Function('u')(x,y)

u=log((x**2+y**2)/(x+y))

xux=simplify(x*diff(u,x))

display(Eq(x*Derivative(f,x),xux))

yuy=simplify(y*diff(u,y))
display(Eq(y*Derivative(f,y),yuy))

LHS= x*Derivative(f,x)+y*Derivative(f,y)

RHS= simplify(xux+yuy)

display(Eq(LHS,RHS))

2a.

from math import *

from sympy import *

r , t = symbols ('r,t')

r1=4*cos ( t )

r2=5*sin( t )

dr1dt = diff( r1 , t )

dr2dt = diff( r2 , t )

tanphi1=r1/dr1dt

tanphi2=r2/dr2dt

theta= solve ( r1-r2 , t )

tanphi1=tanphi1 . subs(t,theta[0])

tanphi2=tanphi2 . subs(t,theta[0])

phi1= atan ( tanphi1 )

phi2= atan ( tanphi2 )

phi = abs( phi1-phi2 )

print ('Angle between curves is %0.4f radians or '%phi,degrees(phi),'degree' )

from numpy import*

from sympy import*

A=[[25,1,2],[1,4,1],[3,10,4]]

print('The given matrix A is')

display(Matrix(A))

eig_val,eig_vec=linalg.eig(A)

print('Eigen values:\n',eig_val)

print('Eigen vectors:\n',eig_vec)

3a.

import numpy as np

import matplotlib.pyplot as plt

# Function to compute the catenary

def catenary(x, a):

return a * np.cosh(x / a)

# Values

a_value = 2

x_values = np.linspace(-5, 5, 100)

# Plotting the catenary

plt.plot(x_values, y_values, label='Catenary: $y = 2 \cosh(x/2)$')

plt.title('Catenary Curve')

plt.xlabel('x-axis')

plt.ylabel('y-axis')

plt.axhline(0, color='black', linewidth=0.5)

plt.axvline(0, color='black', linewidth=0.5)

plt.grid(color='gray', linestyle='--', linewidth=0.5)

plt.legend()

plt.show()

3b.

from sympy import *

r,z,t=symbols('r,z,t')

U=r*cos(z)*sin(t)

V=r*cos(z)*cos(t)

W=r*sin(z)

Ur=simplify(diff(U,r))

Uz=simplify(diff(U,z))

Ut=simplify(diff(U,t))

Vr=simplify(diff(V,r))

Vz=simplify(diff(V,z))

Vt=simplify(diff(V,t))

Wr=simplify(diff(W,r))

Wz=simplify(diff(W,z))
Wt=simplify(diff(W,t))

M=Matrix([[Ur,Uz,Ut],[Vr,Vz,Vt],[Wr,Wz,Wt]])

J=simplify(det(M))

display(Eq(Determinant(M),J))

4a.

from sympy import *

x,y,t= symbols ('x , y , t ')

x=2*cos(t)

y=2*sin(t)

x1=diff(x,t)

y1=diff(y,t)

dydx=(y1/x1)

display(simplify(dydx))

y2=(diff(dydx,t)*(1/x1))

display(simplify(y2))

rho=simplify(((1+dydx**2)**1.5)/y2)

display( rho )

rho1=abs(rho.subs(t,pi/2))

print ('The radius of curvature is ')

display(simplify(rho1))
4b.

from sympy import *

x,y,t= symbols ('x , y , t ')

x=5*cos(t)

y=5*sin(t)

x1=diff(x,t)

y1=diff(y,t)

dydx=(y1/x1)

display(simplify(dydx))

y2=(diff(dydx,t)*(1/x1))

display(simplify(y2))

rho=simplify(((1+dydx**2)**1.5)/y2)

display( rho )

rho1=abs(rho.subs(t,pi/2))

print ('The radius of curvature is ')

display(simplify(rho1))
5a.

from pylab import *

theta= linspace(0,2*pi,1000)

r=2*abs(cos(2*theta))

polar(theta,r)

title("$r=2|cos(2*\\theta))|$")

axvline(color='r')

show()

5b.

from numpy import *

A=[[6,1,2],[1,10,-1],[2,1,-4]]

X=[1,1,1]

lambda_i=0

while(True):

Y=dot(A,X)

lamb=abs(Y).max()

X=Y/lamb

if abs(lamb-lambda_i)<0.001:

break
 else:

lambda_i=lamb

print("\nThe largest Eigen Value is: %.4f"%lambda_i)

print("\n Corresponding Eigen vector is:",X.round(decimals=4))

print("\n The solution is:\n")

print("The largest Eigenvalue is: %.4f"%lamb)

print("Corresponding EigenVector is:",X.round(decimals=4))

6a.

from sympy import *

x= Symbol('x')

y=sin(x)+cos(x)

ys=0

n=int(input("Enter the degree of polynomial required:"))

for i in range(0,n+1):

ys+=x**i*diff(y,x,i).subs({x:0})/factorial(i)

print("Maclaurin's Series is:")

display(Eq(y,ys))

ys=lambdify(x,ys)

from pylab import *

x=linspace(0,pi/2,1000)

plot(x,sin(x)+cos(x),color='red')

plot(x,ys(x),color='green')

axvline(x=0)a

axhline(y=0)
show()

6b.

from sympy import *

A= Matrix([[1,2,-1],[2,1,4],[3,3,4]])

rA=A.rank()

dim=shape(A)

print("The coefficient of the matrix is:")

display(A)

print(f'The rank of the matrix is {rA}')

if rA==dim[1]:

print("System has trivial solution")

else:

print("System has non trivial solution")

7a.

from sympy import *

x,y,z= symbols('x,y,z')

U=x+3*y*y-z*z*z

V=4*x*x*y*z

W=2*z*z-x*y

Ux=simplify(diff(U,x))

Uy=simplify(diff(U,y))

Uz=simplify(diff(U,z))

Vx=simplify(diff(V,x))
Vy=simplify(diff(V,y))

Vz=simplify(diff(V,z))

Wx=simplify(diff(W,x))

Wy=simplify(diff(W,y))

Wz=simplify(diff(W,z))

A=Matrix([[Ux,Uy,Uz],[Vx,Vy,Vz],[Wx,Wy,Wz]])

J=simplify(det(A))

display(Eq(Determinant(A),J))

7b.

from sympy import *

r,t,a,n= symbols ('r , t , a , n')

r=a*sin(n*t)

r1= diff (r , t )

display(r1)

r2= diff ( r1 , t )

display(r2)

rho =(( r**2+r1**2 )**( 1.5 ))/( r**2+2*r1**2-r*r2 ) ;

display(rho)

rho1 = rho . subs (t ,pi/2 )

rho1 = rho . subs (n ,1 )

rho1=simplify(rho1)

print ('The radius of curvature is ' )

display(simplify(rho1))

print('The curvature is')

display(simplify(1/rho1))
8a)

import numpy as np

# Coefficient matrix A

A = np.array([[2, 4, 4],

[4, 2, 2],

[4, 2, 5]])

# Right-hand side vector B

B = np.array([8, 22, 5])

# Use numpy's linear algebra solver to find the solution

try:

solution = np.linalg.solve(A, B)

print("The system has a unique solution:")

print("x =", solution[0])

print("y =", solution[1])

print("z =", solution[2])

except np.linalg.LinAlgError:

print("The system is inconsistent and has no solution.")

8b.

from sympy import *

A=Matrix([[2,1,1],[4,2,-2],[4,22,2]])

B=Matrix([[4],[8],[5]])

dim=shape(A)

AB=A.col_insert(dim[1], B)

rA=A.rank()

rb=AB.rank()

print("The co-efficient matrix A is")

display(A)

print("The augmented matrix AB is")

display(AB)

print(f"rank(A)={rA} and rank(AB)={rb}")

if(rA==rb):

if(rA==dim[1]):

print('System has unique solution')

else:

print('System has infinitely many solution')

print(solve_linear_system(AB,x,y,z))

else:

print('The system of equation is inconsitent')

from sympy import *

x,y,z= symbols('x,y,z')

U=x*y/z

V=y*z/x

W=z*x/y

Ux=simplify(diff(U,x))

Uy=simplify(diff(U,y))

Uz=simplify(diff(U,z))

Vx=simplify(diff(V,x))

Vy=simplify(diff(V,y))

Vz=simplify(diff(V,z))

Wx=simplify(diff(W,x))

Wy=simplify(diff(W,y))

Wz=simplify(diff(W,z))

A=Matrix([[Ux,Uy,Uz],[Vx,Vy,Vz],[Wx,Wy,Wz]])

J=simplify(det(A))

display(Eq(Determinant(A),J))

9b.

import sympy as sp

# Define the variable

t = sp.symbols('t')
# Define the polar curve r(t)

r = -2 * (1 + sp.cos(2 * sp.pi * t))

# Calculate the first and second derivatives of r with respect to t

dr_dt = sp.diff(r, t)

d2r_dt2 = sp.diff(dr_dt, t)

# Calculate the radius of curvature R

R = 1 / sp.Abs(d2r_dt2)

# Calculate the curvature kappa

kappa = 1 / R

# Display the results

print("Radius of Curvature (R):")

display(R)

print("\nCurvature (kappa):")

display(kappa)

10a.

import numpy as np

import matplotlib.pyplot as plt

def cos_maclaurin_series(x, terms=10):

result = 0

for n in range(terms):
 result += (-1)*n * (x*(2*n)) / np.math.factorial(2*n)

return result

# Generate x values in the interval (0, 2*pi)

x_values = np.linspace(0, 2*np.pi, 1000)

# Calculate y values for the Maclaurin series (using 10 terms)

y_maclaurin = cos_maclaurin_series(x_values, terms=10)

# Calculate true values for cos(x)

y_true = np.cos(x_values)

# Plotting

plt.figure(figsize=(8, 6))

plt.plot(x_values, y_true, label='True Function: $y = \cos(x)$', color='blue')

plt.plot(x_values, y_maclaurin, label='Maclaurin Series (10 terms)', linestyle='--', color='red')

plt.title("Maclaurin Series Approximation of $y = \cos(x)$")

plt.xlabel('$x$')

plt.ylabel('$y$')

plt.legend()

plt.grid(True)

plt.show()

10b.

from sympy import *

A= Matrix([[1,2,-1],[2,1,5],[3,3,4]])

rA=A.rank()

dim=shape(A)
print("The coefficient of the matrix is:")

display(A)

print(f'The rank of the matrix is {rA}')

if rA==dim[1]:

print("System has trivial solution")

else:

print("System has non trivial solution")

11a.

import sympy as sp

# Define the variable

t = sp.symbols('t')

# Define the polar curves

r1 = 1 + 2*sp.sin(t)

r2 = 1 + 2*sp.cos(t)

# Define the integrals in the formula

numerator_integral = sp.integrate(r1 * r2, (t, 0, 2*sp.pi))

denominator_integral = sp.sqrt(sp.integrate(r1**2, (t, 0, 2*sp.pi)) * sp.integrate(r2**2, (t, 0,

# Calculate the angle using arccos

angle = sp.acos(numerator_integral / denominator_integral)

# Display the result in degrees

angle_degrees = sp.deg(angle)
print("The angle between the polar curves is approximately:", angle_degrees.evalf(), "degrees.")

11b.

import numpy as np

# Define the matrix A

A = np.array([[11, 1, 2],

[0, 10, 0],

[0, 0, 12]])

# Calculate the eigenvalues and eigenvectors

eigenvalues, eigenvectors = np.linalg.eig(A)

# Display the results

print("Eigenvalues:")

print(eigenvalues)

print("\nEigenvectors:")

print(eigenvectors)

12a.

import sympy as sp

# Define symbols
t, a = sp.symbols('t a')

# Define the polar curves

r1 = sp.exp(a * t)

r2 = sp.exp(-a * t)

# Calculate the angle between the curves

tan_theta = (r2 - r1) / (r1 * r2)

theta = sp.atan(tan_theta)

# Display the result

print("The angle between the polar curves is:")

sp.pprint(sp.deg(theta)) # Convert angle to degrees for better readability

12b.

import numpy as np

# Define the matrix A

A = np.array([[3, 1, 1],

[1, 2, 1],

[1, 1, 12]])

# Calculate eigenvalues and eigenvectors

eigenvalues, eigenvectors = np.linalg.eig(A)

# Display the results

print("Eigenvalues:")
print(eigenvalues)

print("\nEigenvectors:")

print(eigenvectors)

13 a)

import numpy as np

# Coefficient matrix

A = np.array([[1, 1, 1],

[2, 2, -2],

[1, -2, -1]])

# Right-hand side vector

B = np.array([2, 6, 5])

# Solve the system

try:

solution = np.linalg.solve(A, B)

print("Solution:", solution)

except np.linalg.LinAlgError:

print("The system is inconsistent or singular. No unique solution exists.")

13b
import sympy as sp

# Define the symbol x

x = sp.symbols('x')

# Define the expression

expression = sp.log(sp.sin(x)) / (sp.pi/2 - x)**2

# Calculate the limit

limit_result = sp.limit(expression, x, sp.pi/2)

# Display the result

print("Limit as x approaches pi/2:", limit_result)

14a

import numpy as np

def rayleigh_power_iteration(A, X, tol=1e-6, max_iter=1000):

for i in range(max_iter):

Y = np.dot(A, X)

eigenvalue = np.dot(Y.T, X) / np.dot(X.T, X)

eigenvector = Y / np.linalg.norm(Y)

# Check for convergence

if np.linalg.norm(Y - eigenvalue * X) < tol:

return eigenvalue, eigenvector

raise ValueError("Power iteration did not converge within the maximum number of iterations.")

# Given matrix A and initial approximation X

A = np.array([[25, 1, 2],

[1, 3, 0],

[2, 0, -4]])

X = np.array([[1], [0], [1]])

# Compute eigenvalue and eigenvector using Rayleigh power iteration

eigenvalue, eigenvector = rayleigh_power_iteration(A, X)

# Display results

print("Largest Eigenvalue:", eigenvalue)

print("Corresponding Eigenvector:")

print(eigenvector)

14b)

import sympy as sp

# Define symbols

x, y = sp.symbols('x y')

# Define the function u

u = sp.exp(x * sp.cos(y) - y * sp.sin(y))

du2_dx2 = sp.diff(u, x, x)

du2_dy2 = sp.diff(u, y, y)

# Display the results

print("Second partial derivative with respect to x twice (u_xx):", du2_dx2)

print("Second partial derivative with respect to y twice (u_yy):", du2_dy2)

15.a)

import sympy as sp

# Define the variables and the parameter

a, t = sp.symbols('a t')

# Parametric equations

x = a * sp.cos(t)**3

y = a * sp.sin(t)**3

# Calculate the first and second derivatives of y with respect to t

y_prime = sp.diff(y, t)

y_double_prime = sp.diff(y_prime, t)

# Substitute t = pi/4 into the derivatives

t_value = sp.pi / 4

y_prime_value = y_prime.subs(t, t_value)

# Calculate the radius of curvature (R) and curvature (kappa)

R = ((1 + y_prime_value**2)**(3/2)) / sp.Abs(y_double_prime_value)

kappa = sp.Abs(y_double_prime_value) / ((1 + y_prime_value**2)**(3/2))

# Display the results

print(f"At t = π/4:")

print(f"Radius of curvature (R): {R.subs(a, 1)}")

print(f"Curvature (kappa): {kappa.subs(a, 1)}")

15b)

import numpy as np

# Coefficient matrix

A = np.array([[1, 2, -1],

[2, 1, 4],

[1, -1, 5]])

# Calculate the determinant

det_A = np.linalg.det(A)

# Check if the determinant is zero

if det_A == 0:

print("The system has a non-trivial solution.")

else:

print("The system does not have a non-trivial solution.")

import sympy as sp

# Define the symbolic variable

t = sp.symbols('t')

# Given parametric equations

a=1

x = a * (t - sp.sin(t))

y = a * (1 - sp.cos(t))

# Differentiate x and y with respect to t

dx_dt = sp.diff(x, t)

dy_dt = sp.diff(y, t)

# Second derivatives

d2x_dt2 = sp.diff(dx_dt, t)

d2y_dt2 = sp.diff(dy_dt, t)

# Calculate the radius of curvature (R) and curvature (K)

R = ((dx_dt*2 + dy_dt*2)**(3/2)) / abs(dx_dt * d2y_dt2 - d2x_dt2 * dy_dt)

K = abs(dx_dt * d2y_dt2 - d2x_dt2 * dy_dt) / ((dx_dt*2 + dy_dt*2)**(3/2))

# Substitute t=pi/2 into the expressions

t_value = sp.pi / 2

R_value = R.subs(t, t_value)

K_value = K.subs(t, t_value)

print(f"Radius of curvature at t=pi/2: {R_value}")

print(f"Curvature at t=pi/2: {K_value}")

16b)

import numpy as np

# Coefficient matrix

coefficients = np.array([[1, 1, 1],

[0, 1, -1],

[4, -2, -1]])

# Check if the determinant is zero

det = np.linalg.det(coefficients)

if det == 0:

print("The system has non-trivial solutions.")

else:

print("The system only has the trivial solution.")

17a)

import sympy as sp

# Define the variable

# Define the polar curves

r1 = 4 + sp.cos(t)

r2 = -5 * sp.cos(t)

# Define the integrals in the formula

numerator_integral = sp.integrate(r1 * r2, (t, 0, 2*sp.pi))

denominator_integral = sp.sqrt(sp.integrate(r1**2, (t, 0, 2*sp.pi)) * sp.integrate(r2**2, (t, 0,

# Calculate the angle using arccos

angle = sp.acos(numerator_integral / denominator_integral)

# Display the result in degrees

angle_degrees = sp.deg(angle)

print("The angle between the polar curves is approximately:", angle_degrees.evalf(), "degrees.")

17b)

from sympy import *

x,y= symbols('x,y')

f= Function('u')(x,y)

u=log((x**2+y**2)/(x+y))

xux=simplify(x*diff(u,x))

display(Eq(x*Derivative(f,x),xux))

yuy=simplify(y*diff(u,y))

display(Eq(y*Derivative(f,y),yuy))

LHS= x*Derivative(f,x)+y*Derivative(f,y)
RHS= simplify(xux+yuy)

display(Eq(LHS,RHS))

18a)

import numpy as np

import matplotlib.pyplot as plt

# Given parameters

a=4

b=2

# Define the theta values

theta = np.linspace(0, 2 * np.pi, 1000)

# Calculate the corresponding radial values for the ellipse

r = (a * b) / np.sqrt(a*2 * np.sin(theta)**2 + b*2 * np.cos(theta)**2)

# Convert polar coordinates to Cartesian coordinates

x = r * np.cos(theta)

y = r * np.sin(theta)

# Plot the ellipse

plt.figure(figsize=(8, 8))

plt.plot(x, y, label=r'$r = \frac{a \cdot b}{\sqrt{a \cdot \sin^2(\theta) + b \cdot \cos^2(\theta)}}$')

plt.title('Ellipse Plot')

plt.xlabel('x')

plt.ylabel('y')

plt.legend()

plt.grid(True)
plt.axis('equal') # Ensure equal scaling on both axes

plt.show()

18b)

import numpy as np

# Given matrix A

A = np.array([[6, 1, 2],

[1, 10, -1],

[2, 1, -4]])

# Initial approximation to the eigenvector

x = np.array([1, 1, 1])

# Number of iterations

num_iterations = 20

# Rayleigh power method

for i in range(num_iterations):

# Calculate the next approximation of the eigenvector

y = np.dot(A, x)

# Calculate the next approximation of the eigenvalue

eigenvalue = np.dot(x, y) / np.dot(x, x)

# Normalize the eigenvector

x = y / np.linalg.norm(y)
# Display the results

print("Approximated largest eigenvalue:", eigenvalue)

print("Corresponding eigenvector:", x)

19a)

from sympy import *

x= Symbol('x')

y=cos(x)

ys=0

n=int(input("Enter the degree of polynomial required:"))

for i in range(0,n+1):

ys+=x**i*diff(y,x,i).subs({x:0})/factorial(i)

print("Maclaurin's Series is:")

display(Eq(y,ys))

ys=lambdify(x,ys)

from pylab import *

x=linspace(0,pi/2,1000)

plot(x,cos(x),color='red')

plot(x,ys(x),color='green')

axvline(x=0)

axhline(y=0)

show()
19b)

import sympy as sp

# Define the symbolic variable

t = sp.symbols('t')

# Given polar curve

r = 2 * (1 - sp.cos(t))

# Compute the first and second derivatives of r with respect to t

dr_dt = sp.diff(r, t)

d2r_dt2 = sp.diff(dr_dt, t)

# Substitute t = pi/2 into r, dr_dt, and d2r_dt2

t_value = sp.pi / 2

r_value = r.subs(t, t_value)

dr_dt_value = dr_dt.subs(t, t_value)

d2r_dt2_value = d2r_dt2.subs(t, t_value)

# Calculate the radius of curvature (R) and curvature (K)

R = ((r_value*2 + dr_dt_value*2)**(3/2)) / abs(d2r_dt2_value)

K = abs(d2r_dt2_value) / ((r_value*2 + dr_dt_value*2)**(3/2))

# Display the results

print(f"Radius of curvature at t=pi/2: {R.evalf()}")

print(f"Curvature at t=pi/2: {K.evalf()}")

import sympy as sp

# Define symbolic variables

x, y, z = sp.symbols('x y z')

# Given functions u, v, and w

u=x*y/z

v=y*z/x

w=z*x/y

# Compute the partial derivatives

du_dx = sp.diff(u, x)

du_dy = sp.diff(u, y)

du_dz = sp.diff(u, z)

dv_dx = sp.diff(v, x)

dv_dy = sp.diff(v, y)

dv_dz = sp.diff(v, z)

dw_dx = sp.diff(w, x)

dw_dy = sp.diff(w, y)

dw_dz = sp.diff(w, z)

# Create the Jacobian matrix

Jacobian_matrix = sp.Matrix([[du_dx, du_dy, du_dz],

[dv_dx, dv_dy, dv_dz],

# Calculate the Jacobian determinant

Jacobian_det = Jacobian_matrix.det()

# Display the result

print("Jacobian determinant (J):", Jacobian_det.simplify())

20b)

import numpy as np

# Coefficient matrix

A = np.array([[1, 2, -1],

[2, 1, 5],

[3, 3, 4]])

# Calculate the determinant

det_A = np.linalg.det(A)

# Check if the determinant is zero

if det_A == 0:

print("The system has a non-trivial solution.")

else:

print("The system does not have a non-trivial solution.")

import matplotlib.pyplot as plt

# Given points

x = [1, 6, 4, 7, 10]

y = [3, 2, 5, 8, 4]

# Plotting the points

plt.scatter(x, y, color='red', marker='o', label='Data Points')

# Plotting the curve connecting the points

plt.plot(x, y, linestyle='-', color='blue', label='Curve')

# Adding labels and title

plt.xlabel('X-axis')

plt.ylabel('Y-axis')

plt.title('Curve Connecting Points')

# Adding legend

plt.legend()

# Display the plot

plt.show()
21b)

from sympy import *

x,y,z= symbols('x,y,z')

f= Function('u')(x,y,z)

u=log((x**3+y**3+z**3-3*x*y*z))

ux=simplify(diff(u,x))

display(Eq(Derivative(f,x),ux))

uy=simplify(diff(u,y))

display(Eq(Derivative(f,y),uy))

uz=simplify(diff(u,z))

display(Eq(Derivative(f,z),uz))

LHS= Derivative(f,x)+Derivative(f,y)+Derivative(f,z)

RHS= simplify(ux+uy+uz)

display(Eq(LHS,RHS))

22a)

import numpy as np

import matplotlib.pyplot as plt

# Generate theta values

theta = np.linspace(0, 2 * np.pi, 100)

# Parametric equations for a circle

x = 2 * np.cos(theta)

y = 2 * np.sin(theta)
# Plotting the parametric curve (circle)

plt.figure(figsize=(6, 6))

plt.plot(x, y, label='Circle: $x=2\cos(\\theta)$, $y=2\sin(\\theta)$')

plt.title('Parametric Circle')

plt.xlabel('x-axis')

plt.ylabel('y-axis')

plt.axhline(0, color='black',linewidth=0.5)

plt.axvline(0, color='black',linewidth=0.5)

plt.grid(color = 'gray', linestyle = '--', linewidth = 0.5)

plt.legend()

plt.axis('equal') # Equal scaling for x and y axes

plt.show()

22.b)

from sympy import *

r,z,t=symbols('r,z,t')

U=r*cos(z)*sin(t)

V=r*cos(z)*cos(t)

W=r*sin(z)

Ur=simplify(diff(U,r))

Uz=simplify(diff(U,z))

Ut=simplify(diff(U,t))

Vr=simplify(diff(V,r))

Vz=simplify(diff(V,z))

Vt=simplify(diff(V,t))

Wr=simplify(diff(W,r))

Wz=simplify(diff(W,z))
Wt=simplify(diff(W,t))

M=Matrix([[Ur,Uz,Ut],[Vr,Vz,Vt],[Wr,Wz,Wt]])

J=simplify(det(M))

display(Eq(Determinant(M),J))

23a)

import numpy as np

import matplotlib.pyplot as plt

# Function to compute the catenary

def catenary(x, a):

return a * np.cosh(x / a)

# Values

a_value = 2

x_values = np.linspace(-5, 5, 100)

y_values = catenary(x_values, a_value)

# Plotting the catenary

plt.plot(x_values, y_values, label='Catenary: $y = 2 \cosh(x/2)$')

plt.title('Catenary Curve')

plt.xlabel('x-axis')

plt.ylabel('y-axis')

plt.axhline(0, color='black', linewidth=0.5)

plt.axvline(0, color='black', linewidth=0.5)

plt.grid(color='gray', linestyle='--', linewidth=0.5)

plt.legend()
plt.show()

23b)

import sympy as sp

# Define the variables and parameter

t, a = sp.symbols('t a')

a_value = 5

# Define the parametric equations

x = a * sp.cos(t)

y = a * sp.sin(t)

# Calculate the first and second derivatives with respect to t

dx_dt = sp.diff(x, t)

dy_dt = sp.diff(y, t)

d2x_dt2 = sp.diff(dx_dt, t)

d2y_dt2 = sp.diff(dy_dt, t)

# Calculate the curvature formula: |r'(t) x r''(t)| / |r'(t)|^3

curvature = (dx_dt * d2y_dt2 - dy_dt * d2x_dt2) / (dx_dt**2 + dy_dt**2)**(3/2)

# Calculate the radius of curvature formula: 1 / |curvature|

radius_of_curvature = 1 / curvature

# Substitute the given values

curvature_value = curvature.subs({t: t_value, a: a_value})

radius_of_curvature_value = radius_of_curvature.subs({t: t_value, a: a_value})

# Display the results

print("Curvature at t = π/2:", curvature_value.evalf())

print("Radius of Curvature at t = π/2:", radius_of_curvature_value.evalf())

24a)

import numpy as np

def rayleigh_power_iteration(A, X, tol=1e-6, max_iter=1000):

for i in range(max_iter):

Y = np.dot(A, X)

eigenvalue = np.dot(Y.T, X) / np.dot(X.T, X)

eigenvector = Y / np.linalg.norm(Y)

# Check for convergence

if np.linalg.norm(Y - eigenvalue * X) < tol:

return eigenvalue, eigenvector

X = eigenvector

raise ValueError("Power iteration did not converge within the maximum number of iterations.")

# Given matrix A and initial approximation X

X = np.array([[1], [0], [0]])

# Compute eigenvalue and eigenvector using Rayleigh power iteration

eigenvalue, eigenvector = rayleigh_power_iteration(A, X)

# Display results

print("Largest Eigenvalue:", eigenvalue)

print("Corresponding Eigenvector:")

print(eigenvector)

24b)

from sympy import *

x= Symbol('x')

f=((1**(1/x)+2**(1/x)+3**(1/x))/3)**(3*x)

LHS=Limit(f,x,00)

RHS=limit(f,x,00)

display(Eq(LHS,RHS))

25a)

import numpy as np

def is_diagonally_dominant(matrix):
 n = len(matrix)

for i in range(n):

diagonal_sum = sum(abs(matrix[i, j]) for j in range(n) if j != i)

if abs(matrix[i, i]) <= diagonal_sum:

return False

return True

def gauss_seidel(matrix, b, x0=None, tol=1e-6, max_iter=1000):

n = len(matrix)

x = np.zeros(n) if x0 is None else x0.copy()

for _ in range(max_iter):

x_old = x.copy()

for i in range(n):

sigma = np.dot(matrix[i, :i], x[:i]) + np.dot(matrix[i, i+1:], x_old[i+1:])

x[i] = (b[i] - sigma) / matrix[i, i]

# Check for convergence

if np.linalg.norm(x - x_old, np.inf) < tol:

return x

raise ValueError("Gauss-Seidel method did not converge within the maximum number of
iterations.")

# Given system

A = np.array([[10, 1, 1],

[1, 10, 1],

[1, 1, 10]])

b = np.array([12, 12, 12])

if is_diagonally_dominant(A):

# Solve the system using Gauss-Seidel method

solution = gauss_seidel(A, b)

print("Solution using Gauss-Seidel method:", solution)

else:

print("The system is not diagonally dominant. Gauss-Seidel method may not converge.")

26a)

import numpy as np

# Coefficient matrix A

A = np.array([[1,1,1],

[2, 2, -2],

[1, -2, -1]])

# Right-hand side vector B

B = np.array([2,4,5])

# Use numpy's linear algebra solver to find the solution

try:

solution = np.linalg.solve(A, B)

print("The system has a unique solution:")

print("x =", solution[0])

print("z =", solution[2])

except np.linalg.LinAlgError:

print("The system is inconsistent and has no solution.")

26b)

from sympy import *

x= Symbol('x')

f=(1-cos(x))/(x*log(1+x))

LHS=Limit(f,x,00)

RHS=limit(f,x,00)

display(Eq(LHS,RHS))

27a)

import numpy as np

def is_diagonally_dominant(matrix):

n = len(matrix)

for i in range(n):

diagonal_sum = sum(abs(matrix[i, j]) for j in range(n) if j != i)

if abs(matrix[i, i]) <= diagonal_sum:

return False

return True

def gauss_seidel(matrix, b, x0=None, tol=1e-6, max_iter=1000):

n = len(matrix)
 x = np.zeros(n) if x0 is None else x0.copy()

for _ in range(max_iter):

x_old = x.copy()

for i in range(n):

sigma = np.dot(matrix[i, :i], x[:i]) + np.dot(matrix[i, i+1:], x_old[i+1:])

x[i] = (b[i] - sigma) / matrix[i, i]

# Check for convergence

if np.linalg.norm(x - x_old, np.inf) < tol:

return x

raise ValueError("Gauss-Seidel method did not converge within the maximum number of
iterations.")

# Given system

A = np.array([[1,1,1],

[2,2,-2],

[1,-2,-1]])

b = np.array([2,4,5])

# Check for diagonal dominance

if is_diagonally_dominant(A):

# Solve the system using Gauss-Seidel method

solution = gauss_seidel(A, b)

print("Solution using Gauss-Seidel method:", solution)

else:

print("The system is not diagonally dominant. Gauss-Seidel method may not converge.")
28a)

import numpy as np

def is_diagonally_dominant(matrix):

n = len(matrix)

for i in range(n):

diagonal_sum = sum(abs(matrix[i, j]) for j in range(n) if j != i)

if abs(matrix[i, i]) <= diagonal_sum:

return False

return True

def gauss_seidel(matrix, b, x0=None, tol=1e-6, max_iter=1000):

n = len(matrix)

x = np.zeros(n) if x0 is None else x0.copy()

for _ in range(max_iter):

x_old = x.copy()

for i in range(n):

sigma = np.dot(matrix[i, :i], x[:i]) + np.dot(matrix[i, i+1:], x_old[i+1:])

x[i] = (b[i] - sigma) / matrix[i, i]

# Check for convergence

if np.linalg.norm(x - x_old, np.inf) < tol:

return x

raise ValueError("Gauss-Seidel method did not converge within the maximum number of
iterations.")
# Given system

A = np.array([[25,1,1],

[2,10,-3],

[4,-2,-12]])

b = np.array([27,9,-10])

# Check for diagonal dominance

if is_diagonally_dominant(A):

# Solve the system using Gauss-Seidel method

solution = gauss_seidel(A, b)

print("Solution using Gauss-Seidel method:", solution)

else:

print("The system is not diagonally dominant. Gauss-Seidel method may not converge.")

29a)

import numpy as np

# Coefficient matrix A

A = np.array([[25,1,1],

[2,10,-3],

[4,-2,-12]])

# Right-hand side vector B

B = np.array([27,9,-10])

# Use numpy's linear algebra solver to find the solution

try:
 solution = np.linalg.solve(A, B)

print("The system has a unique solution:")

print("x =", solution[0])

print("y =", solution[1])

print("z =", solution[2])

except np.linalg.LinAlgError:

print("The system is inconsistent and has no solution.")

29b)

from sympy import *

x= Symbol('x')

f=((1**(1/x)+2**(1/x)+3**(1/x))/3)**(3*x)

LHS=Limit(f,x,00)

RHS=limit(f,x,00)

display(Eq(LHS,RHS))

30a)

import numpy as np

def rayleigh_power_iteration(A, X, tol=1e-6, max_iter=1000):

for i in range(max_iter):

Y = np.dot(A, X)
 eigenvalue = np.dot(Y.T, X) / np.dot(X.T, X)

eigenvector = Y / np.linalg.norm(Y)

# Check for convergence

if np.linalg.norm(Y - eigenvalue * X) < tol:

return eigenvalue, eigenvector

X = eigenvector

raise ValueError("Power iteration did not converge within the maximum number of iterations.")

# Given matrix A and initial approximation X

A = np.array([[5, 1, 1], [1, 3, -1], [2, -1, -4]])

X = np.array([[1], [0], [0]])

# Compute eigenvalue and eigenvector using Rayleigh power iteration

eigenvalue, eigenvector = rayleigh_power_iteration(A, X)

# Display results

print("Largest Eigenvalue:", eigenvalue)

print("Corresponding Eigenvector:")

print(eigenvector)

30b)

import numpy as np

import matplotlib.pyplot as plt

# Parameter values
a_values = np.linspace(0, 6 * np.pi, 1000)

# Parametric equations

x_values = 2 * (a_values - np.sin(a_values))

y_values = 2 * (1 - np.cos(a_values))

# Plot the parametric curve

plt.plot(x_values, y_values, label='Parametric Curve')

plt.title('Parametric Curve: x = 2(a - sin(a)), y = 2(1 - cos(a))')

plt.xlabel('x')

plt.ylabel('y')

plt.legend()

plt.grid(True)

plt.show()