# Implementation of Source Panel Method as Described in Katz Plotkin Book

I am currently trying to implement Source Panel method as described in Katz and Plotkin in Low Speed Aerodynamics. I have successfully implemented two previous methods. However, I am fully blocked on the source panel method in implementing the SOR2DC function as is called in the book in section 11.2.1.

Theta 1 and Theta 2 Varaible are the same as in the image the same can be said about R1 and R2. P is Point_Collocation.

def Constant_Source_Induced_Vel(Panel_Start_End,Point_Collocation,Source_Strength=1):
Panel_Start=Panel_Start_End[0]
Panel_End=Panel_Start_End[1]
print(Panel_Start)
print(Panel_End)

V_Panel=Panel_End-Panel_Start

V_R1=Point_Collocation-Panel_Start
V_R2=Point_Collocation-Panel_End

R1=np.linalg.norm(V_R1)
R2=np.linalg.norm(V_R2)
print(R1,R2,'Distance')

TH1=np.arccos(np.dot(V_R1,V_Panel)/(np.linalg.norm(V_R1)*np.linalg.norm(V_Panel)))
TH2=np.arccos(np.dot(V_R2,V_Panel)/(np.linalg.norm(V_R2)*np.linalg.norm(V_Panel)))
print(TH1,TH2,'Angles\n')

up=(Source_Strength/(4*np.pi))*np.log(R1**2/R2**2)
wp=(Source_Strength/(2*np.pi))*(TH2-TH1)
return np.array([up,wp])


This is my current implementation I am not using any of the change of reference system transformations to make the troubleshooting easier.Furthermore, it states in the text that this can be done. My issue currently is that the diagonal of the coefficient matrix should be 0.5 but it is not.This is due to the imposition of a boundary condition which make the surface not be crossed by the flow.

My question is is there any way to make this work?(Make the diagonal 0.5)

I have tried changing the way of obtaining the angles using the arctan instead of the dot product. Furthermore, I have tried changing the reference frame but this seems to cause more issues.

Thanks,

PS.

I have edited the code to use the local frame.

def Constant_Source_Induced_Vel(Panel_Start_End, Point_Collocation, Source_Strength=1):
Panel_Start = Panel_Start_End[0]
Panel_End = Panel_Start_End[1]

a=np.arctan2(Panel_End[1]-Panel_Start[1],Panel_End[0]-Panel_Start[0])#Panel Angle

#Local Frame
P_Coll_LxT=Point_Collocation[0]-Panel_Start[0] #Translate x cord Coll Point
P_Coll_LyT=Point_Collocation[1]-Panel_Start[1] #Translate y cord Coll Point

P_Coll_Lx=P_Coll_LxT*np.cos(a) + P_Coll_LyT*np.sin(a)#Rotate to Local Frame
P_Coll_Ly=-P_Coll_LxT*np.sin(a) + P_Coll_LyT*np.cos(a)#Rotate to Local Frame

P_End_LxT=Panel_End[0]-Panel_Start[0] #Translate x cord Panel End Point
P_End_LyT=Panel_End[1]-Panel_Start[1] #Translate y cord Panel End Point

P_End_Lx=P_End_LxT*np.cos(a) +P_End_LyT*np.sin(a)#Rotate to Local Frame
P_End_Ly=0

#Calculate R1,R2,TH1,TH2
R1=(P_Coll_Lx**2 + P_Coll_Ly**2)**0.5
R2=((P_Coll_Lx-P_End_Lx)**2 + P_Coll_Ly**2)**0.5
TH1=np.arctan2(P_Coll_Ly,P_Coll_Lx)
TH2=np.arctan2(P_Coll_Ly,P_Coll_Lx-P_End_Lx)

if Point_Collocation[0] == (Panel_Start[0]+Panel_End[0])/2 and Point_Collocation[1] == (Panel_Start[1]+Panel_End[1])/2:
up = 0
wp = 0.5
print('Diagonal')
print(TH1,TH2,'Angles')
print(R1,R2,'Distances')
else:
up=(Source_Strength/(2*np.pi))*np.log(R1/R2)
wp=(Source_Strength/(2*np.pi))*(TH2-TH1)

u = up *np.cos(a) - wp*np.sin(a)
w = up *np.sin(a) + wp*np.cos(a)
print(u,w,'Final Velocity')
return np.array([u, w])


It results in the following matrix.

Although it seems that the matrix is now ok the issues persists as with an increase in the number of panels the sum of source strenghts does not get closer to zero.

The panels in this graph start from rear bottom and go in a clockwise direction towards the top surface.The pressure on top and bottom should be the same since this is a symmetric airfoil.

• I'm not sure what you hope we can provide in terms of help. Your question is not specific enough to actually contain a question (i.e., there is no question mark), and this isn't the kind of site where you can provide a link to several hundred or thousand lines of Python and hope that someone will go through them to point out your bug. Mar 10 at 11:42
• I think that the panel method is the name given to the boundary element method in potential flow in. That's not my area, though. So, take it with a grain of salt. Mar 10 at 18:29
• Start with figuring out why some of the entries are NaNs. That should be easy enough to track down. Mar 11 at 9:47
• I have added another method which should also work. Mar 11 at 13:18
• This contains the code which I am trying to replicated in the fortan. github.com/alwinw/KatzPlotkinPy/blob/master/tests/f77/SOR2DC.f Mar 11 at 16:56