View raw

1 %wing shear flow 2 clear all; 3 close all; 4 5 6 7 8 Vx = 1; Vz = 1; My = 1; %test loads will be applied individually 9 10 %define a few 11 numTopStringers = 6; 12 numBottomStringers = 8; 13 numNoseTopStringers = 4; 14 numNoseBottomStringers = 4; 15 16 t_upper = 0.02/12; 17 t_lower = 0.02/12; 18 t_upper_front = 0.02/12; 19 t_lower_front = 0.02/12; 20 t_frontSpar = 0.04/12; 21 t_rearSpar = 0.04/12; 22 23 frontSpar = 0.2; 24 backSpar = 0.7; 25 chord = 5; 26 27 sparCaps(1).posX = frontSpar*chord; 28 sparCaps(2).posX = frontSpar*chord; 29 sparCaps(3).posX = backSpar*chord; 30 sparCaps(4).posX = backSpar*chord; 31 32 sparCaps(1).posZ = get_z(frontSpar,1)*chord; 33 sparCaps(2).posZ = get_z(frontSpar,0)*chord; 34 sparCaps(3).posZ = get_z(backSpar,1)*chord; 35 sparCaps(4).posZ = get_z(backSpar,0)*chord; 36 37 sparCaps(1).area = .1; 38 sparCaps(2).area = .1; 39 sparCaps(3).area = .1; 40 sparCaps(4).area = .1; 41 42 upperStringerGap = (sparCaps(3).posX - sparCaps(1).posX)/(numTopStringers + 1); 43 lowerStringerGap = (sparCaps(3).posX - sparCaps(1).posX)/(numBottomStringers + 1); 44 upperNoseStringerGap = (sparCaps(1).posX - 0)/(numNoseTopStringers + 1); 45 lowerNoseStringerGap = (sparCaps(1).posX - 0)/(numNoseBottomStringers + 1); 46 47 48 %set stringers spaced evenly along X axis betwen Spars 49 %top Stringers 50 for i=1:numTopStringers 51 topStringers(i).posX = sparCaps(1).posX + upperStringerGap*i; 52 topStringers(i).posZ = get_z(topStringers(i).posX/chord,1)*chord; 53 topStringers(i).area = .1; 54 end 55 56 %bottom Stringers 57 for i=1:numBottomStringers 58 bottomStringers(i).posX = sparCaps(4).posX - lowerStringerGap*i; 59 bottomStringers(i).posZ = get_z(bottomStringers(i).posX/chord,0)*chord; 60 bottomStringers(i).area = .1; 61 62 end 63 64 %nose bottom Stringers 65 for i=1:numNoseBottomStringers 66 noseBottomStringers(i).posX = sparCaps(2).posX - lowerNoseStringerGap*i; 67 noseBottomStringers(i).posZ = get_z(noseBottomStringers(i).posX/chord,0)*chord; 68 noseBottomStringers(i).area = .1; 69 end 70 71 %nose top Stringers 72 for i=1:numNoseTopStringers 73 noseTopStringers(i).posX = upperNoseStringerGap*i; 74 noseTopStringers(i).posZ = get_z(noseTopStringers(i).posX/chord,1)*chord; 75 noseTopStringers(i).area = .1; 76 end 77 78 79 centroid.posX = sum([sparCaps.posX].*[sparCaps.area]) + ... 80 sum([topStringers.posX].*[topStringers.area]) + ... 81 sum([bottomStringers.posX].*[bottomStringers.area]) + ... 82 sum([noseTopStringers.posX].*[noseTopStringers.area]) + ... 83 sum([noseBottomStringers.posX].*[noseBottomStringers.area]); 84 85 centroid.posX = centroid.posX / ( sum([sparCaps.area]) + sum([topStringers.area]) + ... 86 sum([bottomStringers.area]) + sum([noseTopStringers.area]) + sum([noseBottomStringers.area])); 87 88 centroid.posZ = sum([sparCaps.posZ].*[sparCaps.area]) + ... 89 sum([topStringers.posZ].*[topStringers.area]) + ... 90 sum([bottomStringers.posZ].*[bottomStringers.area]) + ... 91 sum([noseTopStringers.posZ].*[noseTopStringers.area]) + ... 92 sum([noseBottomStringers.posZ].*[noseBottomStringers.area]); 93 94 centroid.posZ = centroid.posZ / ( sum([sparCaps.area]) + sum([topStringers.area]) + ... 95 sum([bottomStringers.area]) + sum([noseTopStringers.area]) + sum([noseBottomStringers.area])); 96 97 %summing contributions for inertia terms 98 Ix = 0; Iz = 0; Ixz = 0; 99 100 for i=1:4 %spar caps 101 Ix = Ix + sparCaps(i).area*(sparCaps(i).posZ-centroid.posZ)^2; 102 Iz = Iz + sparCaps(i).area*(sparCaps(i).posX-centroid.posX)^2; 103 Ixz = Ixz + sparCaps(i).area*(sparCaps(i).posX-centroid.posX)*(sparCaps(i).posZ-centroid.posZ); 104 end 105 106 107 for i=1:numTopStringers %top stringers 108 Ix = Ix + topStringers(i).area*(topStringers(i).posZ-centroid.posZ)^2; 109 Iz = Iz + topStringers(i).area*(topStringers(i).posX-centroid.posX)^2; 110 Ixz = Ixz + topStringers(i).area*(topStringers(i).posX-centroid.posX)*(topStringers(i).posZ-centroid.posZ); 111 end 112 for i=1:numBottomStringers %bottom stringers 113 Ix = Ix + bottomStringers(i).area*(bottomStringers(i).posZ-centroid.posZ)^2; 114 Iz = Iz + bottomStringers(i).area*(bottomStringers(i).posX-centroid.posX)^2; 115 Ixz = Ixz + bottomStringers(i).area*(bottomStringers(i).posX-centroid.posX)*(bottomStringers(i).posZ-centroid.posZ); 116 end 117 for i=1:numNoseTopStringers %nose top stringers 118 Ix = Ix + noseTopStringers(i).area*(noseTopStringers(i).posZ-centroid.posZ)^2; 119 Iz = Iz + noseTopStringers(i).area*(noseTopStringers(i).posX-centroid.posX)^2; 120 Ixz = Ixz + noseTopStringers(i).area*(noseTopStringers(i).posX-centroid.posX)*(noseTopStringers(i).posZ-centroid.posZ); 121 end 122 for i=1:numNoseBottomStringers %nose bottom stringers 123 Ix = Ix + noseBottomStringers(i).area*(noseBottomStringers(i).posZ-centroid.posZ)^2; 124 Iz = Iz + noseBottomStringers(i).area*(noseBottomStringers(i).posX-centroid.posX)^2; 125 Ixz = Ixz + noseBottomStringers(i).area*(noseBottomStringers(i).posX-centroid.posX)*(noseBottomStringers(i).posZ-centroid.posZ); 126 end 127 128 %Ixz = -Ixz; 129 130 %define webs 131 132 %% web cell 1 133 134 %upper webs 135 numStringers = numTopStringers; 136 stringerGap = upperStringerGap; 137 webThickness = t_upper; 138 tempStringers = topStringers; 139 140 for i=1:(numStringers+1) 141 web(i).xStart = sparCaps(1).posX + stringerGap*(i-1); 142 web(i).xEnd = sparCaps(1).posX + stringerGap*(i); 143 web(i).thickness = webThickness; 144 web(i).zStart = get_z(web(i).xStart/chord,1)*chord; 145 web(i).zEnd = get_z(web(i).xEnd/chord,1)*chord; 146 if i==1 147 web(i).dp_area = sparCaps(1).area; 148 web(i).dP_X = 0; 149 web(i).dP_Z = 0; 150 web(i).qPrime_X = 0; 151 web(i).qPrime_Z = 0; 152 else 153 web(i).dp_area = tempStringers(i-1).area; 154 dx = web(i).xStart-centroid.posX; dz = web(i).zStart-centroid.posZ; 155 web(i).dP_X = get_dp(dx,dz,Vx,0,Ix,Iz,Ixz,web(i).dp_area); %just Vx 156 web(i).dP_Z = get_dp(dx,dz,0,Vz,Ix,Iz,Ixz,web(i).dp_area); %just Vz 157 web(i).qPrime_X = web(i-1).qPrime_X - web(i).dP_X; 158 web(i).qPrime_Z = web(i-1).qPrime_Z - web(i).dP_Z; 159 end 160 tempInt = get_int(web(i).xStart/chord,web(i).xEnd/chord,1)*chord^2; %integral of airfoil function 161 triangle1 = abs( (web(i).xStart - sparCaps(1).posX)*web(i).zStart/2); 162 triangle2 = abs((web(i).xEnd - sparCaps(1).posX)*web(i).zEnd/2); 163 web(i).Area = tempInt + triangle1 - triangle2; 164 web(i).ds = get_ds(web(i).xStart/chord,web(i).xEnd/chord,1)*chord; 165 web(i).dS_over_t = web(i).ds / web(i).thickness; 166 167 web(i).q_dS_over_t_X = web(i).qPrime_X * web(i).dS_over_t; 168 web(i).q_dS_over_t_Z = web(i).qPrime_Z * web(i).dS_over_t; 169 web(i).two_A_qprime_X = 2*web(i).Area*web(i).qPrime_X; 170 web(i).two_A_qprime_Z = 2*web(i).Area*web(i).qPrime_Z; 171 web(i).qp_dx_X = web(i).qPrime_X *(web(i).xEnd-web(i).xStart); 172 web(i).qp_dx_Z = web(i).qPrime_Z *(web(i).xEnd-web(i).xStart); 173 web(i).qp_dz_X = web(i).qPrime_X *(web(i).zEnd-web(i).zStart); 174 web(i).qp_dz_Z = web(i).qPrime_Z *(web(i).zEnd-web(i).zStart); 175 end 176 webTop = web; 177 web = []; 178 179 %rear spar 180 i=1; 181 web(i).xStart = sparCaps(3).posX; 182 web(i).xEnd = sparCaps(4).posX; 183 web(i).thickness = t_rearSpar; 184 web(i).zStart = sparCaps(3).posZ; 185 web(i).zEnd = sparCaps(4).posZ; 186 web(i).dp_area = sparCaps(3).area; 187 dx = web(i).xStart-centroid.posX; dz = web(i).zStart-centroid.posZ; 188 web(i).dP_X = get_dp(dx,dz,Vx,0,Ix,Iz,Ixz,web(i).dp_area); 189 web(i).dP_Z = get_dp(dx,dz,0,Vz,Ix,Iz,Ixz,web(i).dp_area); 190 web(i).qPrime_X = webTop(numTopStringers+1).qPrime_X - web(i).dP_X; 191 web(i).qPrime_Z = webTop(numTopStringers+1).qPrime_Z - web(i).dP_Z; 192 193 web(i).Area = (sparCaps(3).posX-sparCaps(1).posX)*sparCaps(3).posZ/2 + ... 194 abs((sparCaps(3).posX-sparCaps(1).posX)*sparCaps(4).posZ/2); 195 web(i).ds = abs(sparCaps(3).posZ - sparCaps(4).posZ); 196 web(i).dS_over_t = web(i).ds / web(i).thickness; 197 198 web(i).q_dS_over_t_X = web(i).qPrime_X * web(i).dS_over_t; 199 web(i).q_dS_over_t_Z = web(i).qPrime_Z * web(i).dS_over_t; 200 web(i).two_A_qprime_X = 2*web(i).Area*web(i).qPrime_X; 201 web(i).two_A_qprime_Z = 2*web(i).Area*web(i).qPrime_Z; 202 web(i).qp_dx_X = web(i).qPrime_X *(web(i).xEnd-web(i).xStart); 203 web(i).qp_dx_Z = web(i).qPrime_Z *(web(i).xEnd-web(i).xStart); 204 web(i).qp_dz_X = web(i).qPrime_X *(web(i).zEnd-web(i).zStart); 205 web(i).qp_dz_Z = web(i).qPrime_Z *(web(i).zEnd-web(i).zStart); 206 207 webRearSpar = web; 208 web = []; 209 210 211 %lower webs 212 numStringers = numBottomStringers; 213 stringerGap = lowerStringerGap; 214 webThickness = t_lower; 215 tempStringers = bottomStringers; 216 217 for i=1:(numStringers+1) 218 web(i).xStart = sparCaps(4).posX - stringerGap*(i-1); 219 web(i).xEnd = sparCaps(4).posX - stringerGap*(i); 220 web(i).thickness = webThickness; 221 web(i).zStart = get_z(web(i).xStart/chord,0)*chord; 222 web(i).zEnd = get_z(web(i).xEnd/chord,0)*chord; 223 dx = web(i).xStart-centroid.posX; dz = web(i).zStart-centroid.posZ; 224 if i==1 225 web(i).dp_area = sparCaps(4).area; 226 web(i).dP_X = get_dp(dx,dz,Vx,0,Ix,Iz,Ixz,web(i).dp_area); 227 web(i).dP_Z = get_dp(dx,dz,0,Vz,Ix,Iz,Ixz,web(i).dp_area); 228 web(i).qPrime_X = webRearSpar.qPrime_X - web(i).dP_X; 229 web(i).qPrime_Z = webRearSpar.qPrime_Z - web(i).dP_Z; 230 else 231 web(i).dp_area = tempStringers(i-1).area; 232 web(i).dP_X = get_dp(dx,dz, Vx,0,Ix,Iz,Ixz,web(i).dp_area); 233 web(i).dP_Z = get_dp(dx,dz, 0,Vz,Ix,Iz,Ixz,web(i).dp_area); 234 web(i).qPrime_X = web(i-1).qPrime_X - web(i).dP_X; 235 web(i).qPrime_Z = web(i-1).qPrime_Z - web(i).dP_Z; 236 end 237 238 tempInt = get_int(web(i).xEnd/chord,web(i).xStart/chord,0)*chord^2; %integral of airfoil function 239 triangle2 = abs((web(i).xStart - sparCaps(1).posX)*web(i).zStart/2); 240 triangle1 = abs((web(i).xEnd - sparCaps(1).posX)*web(i).zEnd/2); 241 web(i).Area = tempInt + triangle1 - triangle2; 242 web(i).ds = get_ds(web(i).xStart/chord,web(i).xEnd/chord,0)*chord; 243 web(i).dS_over_t = web(i).ds / web(i).thickness; 244 245 web(i).q_dS_over_t_X = web(i).qPrime_X * web(i).dS_over_t; 246 web(i).q_dS_over_t_Z = web(i).qPrime_Z * web(i).dS_over_t; 247 web(i).two_A_qprime_X = 2*web(i).Area*web(i).qPrime_X; 248 web(i).two_A_qprime_Z = 2*web(i).Area*web(i).qPrime_Z; 249 web(i).qp_dx_X = web(i).qPrime_X*(web(i).xEnd-web(i).xStart); 250 web(i).qp_dx_Z = web(i).qPrime_Z*(web(i).xEnd-web(i).xStart); 251 web(i).qp_dz_X = web(i).qPrime_X*(web(i).zEnd-web(i).zStart); 252 web(i).qp_dz_Z = web(i).qPrime_Z*(web(i).zEnd-web(i).zStart); 253 254 %web(i).radCurv = ... Example: get_curve(web(i).xStart,web(i).xEnd,1) 255 end 256 webBottom = web; 257 web = []; 258 259 %front Spar 260 i=1; 261 web(i).xStart = sparCaps(2).posX; 262 web(i).xEnd = sparCaps(1).posX; 263 web(i).thickness = t_frontSpar; 264 web(i).zStart = sparCaps(2).posZ; 265 web(i).zEnd = sparCaps(1).posZ; 266 web(i).dp_area = sparCaps(2).area; 267 dx = web(i).xStart-centroid.posX; dz = web(i).zStart-centroid.posZ; 268 web(i).dP_X = get_dp(dx,dz,Vx,0,Ix,Iz,Ixz,web(i).dp_area); 269 web(i).dP_Z = get_dp(dx,dz,0,Vz,Ix,Iz,Ixz,web(i).dp_area); 270 web(i).qPrime_X = webBottom(numBottomStringers+1).qPrime_X - web(i).dP_X; 271 web(i).qPrime_Z = webBottom(numBottomStringers+1).qPrime_Z - web(i).dP_Z; 272 web(i).Area = 0; 273 web(i).ds = abs(sparCaps(2).posZ - sparCaps(1).posZ); 274 web(i).dS_over_t = web(i).ds / web(i).thickness; 275 276 web(i).q_dS_over_t_X = web(i).qPrime_X * web(i).dS_over_t; 277 web(i).q_dS_over_t_Z = web(i).qPrime_Z * web(i).dS_over_t; 278 web(i).two_A_qprime_X = 2*web(i).Area*web(i).qPrime_X; 279 web(i).two_A_qprime_Z = 2*web(i).Area*web(i).qPrime_Z; 280 web(i).qp_dx_X = web(i).qPrime_X *(web(i).xEnd-web(i).xStart); 281 web(i).qp_dx_Z = web(i).qPrime_Z *(web(i).xEnd-web(i).xStart); 282 web(i).qp_dz_X = web(i).qPrime_X *(web(i).zEnd-web(i).zStart); 283 web(i).qp_dz_Z = web(i).qPrime_Z *(web(i).zEnd-web(i).zStart); 284 285 webFrontSpar = web; 286 web = []; 287 288 289 290 291 %% web cell 2 292 293 %lower nose webs 294 numStringers = numNoseBottomStringers; 295 stringerGap = lowerNoseStringerGap; 296 webThickness = t_lower_front; 297 tempStringers = noseBottomStringers; 298 299 for i=1:(numStringers+1) 300 web(i).xStart = sparCaps(2).posX - stringerGap*(i-1); 301 web(i).xEnd = sparCaps(2).posX - stringerGap*(i); 302 web(i).thickness = webThickness; 303 web(i).zStart = get_z(web(i).xStart/chord,0)*chord; 304 web(i).zEnd = get_z(web(i).xEnd/chord,0)*chord; 305 dx = web(i).xStart-centroid.posX; dz = web(i).zStart-centroid.posZ; 306 307 if i==1 308 web(i).dp_area = sparCaps(2).area; 309 web(i).dP_X = 0; 310 web(i).dP_Z = 0; 311 web(i).qPrime_X = 0; 312 web(i).qPrime_Z = 0; 313 else 314 web(i).dp_area = tempStringers(i-1).area; 315 web(i).dP_X = get_dp(dx,dz,Vx,0,Ix,Iz,Ixz,web(i).dp_area); 316 web(i).dP_Z = get_dp(dx,dz,0,Vz,Ix,Iz,Ixz,web(i).dp_area); 317 web(i).qPrime_X = web(i-1).qPrime_X - web(i).dP_X; 318 web(i).qPrime_Z = web(i-1).qPrime_Z - web(i).dP_Z; 319 end 320 tempInt = get_int(web(i).xEnd/chord,web(i).xStart/chord,0)*chord^2; %integral of airfoil function 321 triangle1 = abs((web(i).xStart - sparCaps(2).posX)*web(i).zStart/2); 322 triangle2 = abs((web(i).xEnd - sparCaps(2).posX)*web(i).zEnd/2); 323 web(i).Area = tempInt + triangle1 - triangle2; 324 web(i).ds = get_ds(web(i).xStart/chord,web(i).xEnd/chord,0)*chord; 325 web(i).dS_over_t = web(i).ds / web(i).thickness; 326 327 web(i).q_dS_over_t_X = web(i).qPrime_X * web(i).dS_over_t; 328 web(i).q_dS_over_t_Z = web(i).qPrime_Z * web(i).dS_over_t; 329 web(i).two_A_qprime_X = 2*web(i).Area*web(i).qPrime_X; 330 web(i).two_A_qprime_Z = 2*web(i).Area*web(i).qPrime_Z; 331 web(i).qp_dx_X = web(i).qPrime_X *(web(i).xEnd-web(i).xStart); 332 web(i).qp_dx_Z = web(i).qPrime_Z *(web(i).xEnd-web(i).xStart); 333 web(i).qp_dz_X = web(i).qPrime_X *(web(i).zEnd-web(i).zStart); 334 web(i).qp_dz_Z = web(i).qPrime_Z *(web(i).zEnd-web(i).zStart); 335 336 %web(i).radCurv = ... Example: get_curve(web(i).xStart,web(i).xEnd,1) 337 end 338 webLowerNose = web; 339 web = []; 340 341 %upper nose webs 342 numStringers = numNoseTopStringers; 343 stringerGap = upperNoseStringerGap; 344 webThickness = t_upper_front; 345 tempStringers = noseTopStringers; 346 347 for i=1:(numStringers+1) 348 web(i).xStart = stringerGap*(i-1); 349 web(i).xEnd = stringerGap*(i); 350 web(i).thickness = webThickness; 351 web(i).zStart = get_z(web(i).xStart/chord,1)*chord; 352 web(i).zEnd = get_z(web(i).xEnd/chord,1)*chord; 353 dx = web(i).xStart-centroid.posX; dz = web(i).zStart-centroid.posZ; 354 if i==1 355 web(i).dp_area = 0; 356 web(i).dP_X = 0; 357 web(i).dP_Z = 0; 358 web(i).qPrime_X = webLowerNose(numNoseBottomStringers+1).qPrime_X - web(i).dP_X; 359 web(i).qPrime_Z = webLowerNose(numNoseBottomStringers+1).qPrime_Z - web(i).dP_Z; 360 else 361 web(i).dp_area = tempStringers(i-1).area; 362 web(i).dP_X = get_dp(dx,dz,Vx,0,Ix,Iz,Ixz,web(i).dp_area); 363 web(i).dP_Z = get_dp(dx,dz,0,Vz,Ix,Iz,Ixz,web(i).dp_area); 364 web(i).qPrime_X = web(i-1).qPrime_X - web(i).dP_X; 365 web(i).qPrime_Z = web(i-1).qPrime_Z - web(i).dP_Z; 366 end 367 tempInt = get_int(web(i).xStart/chord,web(i).xEnd/chord,1)*chord^2; %integral of airfoil function 368 triangle2 = abs((web(i).xStart - sparCaps(2).posX)*web(i).zStart/2); 369 triangle1 = abs((web(i).xEnd - sparCaps(2).posX)*web(i).zEnd/2); 370 web(i).Area = tempInt + triangle1 - triangle2; 371 web(i).ds = get_ds(web(i).xStart/chord,web(i).xEnd/chord,1)*chord; 372 web(i).dS_over_t = web(i).ds / web(i).thickness; 373 374 web(i).q_dS_over_t_X = web(i).qPrime_X * web(i).dS_over_t; 375 web(i).q_dS_over_t_Z = web(i).qPrime_Z * web(i).dS_over_t; 376 web(i).two_A_qprime_X = 2*web(i).Area*web(i).qPrime_X; 377 web(i).two_A_qprime_Z = 2*web(i).Area*web(i).qPrime_Z; 378 web(i).qp_dx_X = web(i).qPrime_X *(web(i).xEnd-web(i).xStart); 379 web(i).qp_dx_Z = web(i).qPrime_Z *(web(i).xEnd-web(i).xStart); 380 web(i).qp_dz_X = web(i).qPrime_X *(web(i).zEnd-web(i).zStart); 381 web(i).qp_dz_Z = web(i).qPrime_Z *(web(i).zEnd-web(i).zStart); 382 383 end 384 webUpperNose = web; 385 web = []; 386 387 388 %front Spar 389 i=1; 390 web(i).xStart = sparCaps(1).posX; 391 web(i).xEnd = sparCaps(2).posX; 392 web(i).thickness = t_frontSpar; 393 web(i).zStart = sparCaps(1).posZ; 394 web(i).zEnd = sparCaps(2).posZ; 395 web(i).dp_area = sparCaps(1).area; 396 dx = web(i).xStart-centroid.posX; dz = web(i).zStart-centroid.posZ; 397 398 web(i).dP_X = get_dp(dx,dz,Vx,0,Ix,Iz,Ixz,web(i).dp_area); 399 web(i).dP_Z = get_dp(dx,dz,0,Vz,Ix,Iz,Ixz,web(i).dp_area); 400 web(i).qPrime_X = webUpperNose(numNoseTopStringers+1).qPrime_X - web(i).dP_X; 401 web(i).qPrime_Z = webUpperNose(numNoseTopStringers+1).qPrime_Z - web(i).dP_Z; 402 web(i).Area = 0; 403 web(i).ds = abs(sparCaps(1).posZ - sparCaps(2).posZ); 404 web(i).dS_over_t = web(i).ds / web(i).thickness; 405 web(i).q_dS_over_t_X = web(i).qPrime_X * web(i).dS_over_t; 406 web(i).q_dS_over_t_Z = web(i).qPrime_Z * web(i).dS_over_t; 407 web(i).two_A_qprime_X = 2*web(i).Area*web(i).qPrime_X; 408 web(i).two_A_qprime_Z = 2*web(i).Area*web(i).qPrime_Z; 409 web(i).qp_dx_X = web(i).qPrime_X *(web(i).xEnd-web(i).xStart); 410 web(i).qp_dx_Z = web(i).qPrime_Z *(web(i).xEnd-web(i).xStart); 411 web(i).qp_dz_X = web(i).qPrime_X *(web(i).zEnd-web(i).zStart); 412 web(i).qp_dz_Z = web(i).qPrime_Z *(web(i).zEnd-web(i).zStart); 413 414 webFrontSparCell2 = web; 415 web = []; 416 417 418 %check that q'*dx sums up to Vx 419 420 Fx = sum([webTop.qp_dx_X])+webRearSpar.qp_dx_X+ sum([webBottom.qp_dx_X])+webFrontSpar.qp_dx_X; %cell 1 421 Fx = Fx + sum([webLowerNose.qp_dx_X])+ sum([webUpperNose.qp_dx_X]); %cell 2 422 Fx 423 Fz = sum([webTop.qp_dz_X])+webRearSpar.qp_dz_X+ sum([webBottom.qp_dz_X])+webFrontSpar.qp_dz_X; %cell 1 424 Fz = Fz + sum([webLowerNose.qp_dz_X])+ sum([webUpperNose.qp_dz_X]); %cell 2 425 Fz 426 427 %check that q'*dz sums up to Vz 428 429 430 Fx = sum([webTop.qp_dx_Z])+webRearSpar.qp_dx_Z+ sum([webBottom.qp_dx_Z])+webFrontSpar.qp_dx_Z; %cell 1 431 Fx = Fx + sum([webLowerNose.qp_dx_Z])+ sum([webUpperNose.qp_dx_Z]); %cell 2 432 Fx 433 Fz = sum([webTop.qp_dz_Z])+webRearSpar.qp_dz_Z+ sum([webBottom.qp_dz_Z])+webFrontSpar.qp_dz_Z; %cell 1 434 Fz = Fz + sum([webLowerNose.qp_dz_Z])+ sum([webUpperNose.qp_dz_Z]); %cell 2 435 Fz 436 437 %% 438 439 % sum up the ds/t and q*ds/t to solve 2 equations, 2 unknowns 440 441 % [A]*[q1s q2s] = B 442 443 A11 = sum([webTop.dS_over_t])+webRearSpar.dS_over_t+ sum([webBottom.dS_over_t])+webFrontSpar.dS_over_t; 444 A22 = sum([webLowerNose.dS_over_t])+ sum([webUpperNose.dS_over_t])+webFrontSparCell2.dS_over_t; 445 A12 = -webFrontSpar.dS_over_t; 446 A21 = -webFrontSparCell2.dS_over_t; 447 448 B1_X = sum([webTop.q_dS_over_t_X])+webRearSpar.q_dS_over_t_X+ sum([webBottom.q_dS_over_t_X])+webFrontSpar.q_dS_over_t_X; 449 B2_X = sum([webLowerNose.q_dS_over_t_X])+ sum([webUpperNose.q_dS_over_t_X])+webFrontSparCell2.q_dS_over_t_X; 450 B1_Z = sum([webTop.q_dS_over_t_Z])+webRearSpar.q_dS_over_t_Z+ sum([webBottom.q_dS_over_t_Z])+webFrontSpar.q_dS_over_t_Z; 451 B2_Z = sum([webLowerNose.q_dS_over_t_Z])+ sum([webUpperNose.q_dS_over_t_Z])+webFrontSparCell2.q_dS_over_t_Z; 452 453 Amat = [A11 A12; A21 A22]; 454 Bmat_X = -[B1_X;B2_X]; 455 Bmat_Z = -[B1_Z;B2_Z]; 456 457 qs_X = inv(Amat)*Bmat_X; 458 qs_Z = inv(Amat)*Bmat_Z; 459 460 461 462 sum_2_a_q_X = sum([webTop.two_A_qprime_X])+webRearSpar.two_A_qprime_X+ sum([webBottom.two_A_qprime_X]); %cell 1 qprimes 463 sum_2_a_q_X = sum_2_a_q_X + sum([webLowerNose.two_A_qprime_X])+ sum([webUpperNose.two_A_qprime_X]); %cell 2 qprimes 464 sum_2_a_q_X = sum_2_a_q_X + 2*qs_X(1)*(sum([webTop.Area])+webRearSpar.Area+ sum([webBottom.Area])); 465 sum_2_a_q_X = sum_2_a_q_X + 2*qs_X(2)*(sum([webLowerNose.Area])+ sum([webUpperNose.Area])); 466 467 sum_2_a_q_Z = sum([webTop.two_A_qprime_Z])+webRearSpar.two_A_qprime_Z+ sum([webBottom.two_A_qprime_Z]); %cell 1 qprimes 468 sum_2_a_q_Z = sum_2_a_q_Z + sum([webLowerNose.two_A_qprime_Z])+ sum([webUpperNose.two_A_qprime_Z]); %cell 2 qprimes 469 sum_2_a_q_Z = sum_2_a_q_Z + 2*qs_Z(1)*(sum([webTop.Area])+webRearSpar.Area+ sum([webBottom.Area])); 470 sum_2_a_q_Z = sum_2_a_q_Z + 2*qs_Z(2)*(sum([webLowerNose.Area])+ sum([webUpperNose.Area])); 471 472 %shear center 473 sc.posX = sum_2_a_q_Z / Vz + frontSpar*chord; 474 sc.posZ = - sum_2_a_q_X / Vx; 475 476 477 % now consider the torque representing shifting the load from the quarter 478 % chord to the SC (need to check signs on these moments) 479 480 torque_Z = Vz*(sc.posX - 0.25*chord); 481 torque_X = -Vx*sc.posZ; 482 483 484 Area1 = sum([webTop.Area]) + webRearSpar.Area + sum([webBottom.Area]); 485 %check area 486 Area1_check = get_int(frontSpar,backSpar,1)*chord^2 + get_int(frontSpar,backSpar,0)*chord^2; 487 488 Area2 = sum([webLowerNose.Area]) + sum([webUpperNose.Area]); 489 Area2_check = get_int(0,frontSpar,1)*chord^2 + get_int(0,frontSpar,0)*chord^2; 490 491 492 %for twist equation (see excel spreadsheet example) 493 494 q1t_over_q2t = (A22/Area2 + webFrontSpar.dS_over_t/Area1)/(A11/Area1 + webFrontSpar.dS_over_t/Area2); 495 496 q2t = torque_X/(2*Area1*q1t_over_q2t + 2*Area2); 497 q1t = q2t*q1t_over_q2t; 498 qt_X = [q1t;q2t]; 499 500 q2t = torque_Z/(2*Area1*q1t_over_q2t + 2*Area2); 501 q1t = q2t*q1t_over_q2t; 502 qt_Z = [q1t;q2t]; 503 504 505 506 % --- - add up all shear flows: qtot = (qPrime + qs) + qt 507 508 509 510 511 %--- insert force balance to check total shear flows --- 512 513 % --- -- 514 515 516 %end 517 518 sc 519 520 521 %plotting airfoil cross-section 522 523 xChord = 0:.01:1; 524 xChord = xChord*chord; 525 upperSurface = zeros(1,length(xChord)); 526 lowerSurface = zeros(1,length(xChord)); 527 528 for i=1:length(xChord) 529 upperSurface(i) = get_z(xChord(i)/chord,1)*chord; 530 lowerSurface(i) = get_z(xChord(i)/chord,0)*chord; 531 end 532 533 figure; hold on; axis equal; grid on; 534 %plot(xChord,z_camber,'-') 535 plot(xChord,upperSurface,'-k','linewidth',2) 536 plot(xChord,lowerSurface,'-k','linewidth',2) 537 plot([0 1],[0 0],'--k','linewidth',1) 538 539 540 for i = 1:length(webTop) 541 vecX = [frontSpar*chord webTop(i).xStart webTop(i).xEnd]; 542 vecZ = [0 webTop(i).zStart webTop(i).zEnd]; 543 fill(vecX,vecZ,[0.9 0.9 0.9]) 544 end 545 546 for i = 1:length(webBottom) 547 vecX = [frontSpar*chord webBottom(i).xStart webBottom(i).xEnd]; 548 vecZ = [0 webBottom(i).zStart webBottom(i).zEnd]; 549 fill(vecX,vecZ,[0.9 0.9 0.9]) 550 end 551 552 for i = 1:length(webUpperNose) 553 vecX = [frontSpar*chord webUpperNose(i).xStart webUpperNose(i).xEnd]; 554 vecZ = [0 webUpperNose(i).zStart webUpperNose(i).zEnd]; 555 fill(vecX,vecZ,[0.7 0.9 1.0]) 556 end 557 558 for i = 1:length(webLowerNose) 559 vecX = [frontSpar*chord webLowerNose(i).xStart webLowerNose(i).xEnd]; 560 vecZ = [0 webLowerNose(i).zStart webLowerNose(i).zEnd]; 561 fill(vecX,vecZ,[0.7 0.9 1.0]) 562 end 563 564 vecX = [frontSpar*chord sparCaps(3).posX sparCaps(4).posX]; 565 vecZ = [0 sparCaps(3).posZ sparCaps(4).posZ]; 566 fill(vecX,vecZ,[0.9 0.9 0.9]) 567 568 569 sparCapSize = 18; 570 stringerSize = 18; 571 plot([sparCaps(1).posX sparCaps(2).posX],[sparCaps(1).posZ sparCaps(2).posZ],'-k','linewidth',2) 572 plot([sparCaps(3).posX sparCaps(4).posX],[sparCaps(3).posZ sparCaps(4).posZ],'-k','linewidth',2) 573 plot([sparCaps.posX],[sparCaps.posZ],'.b','markersize',sparCapSize) 574 plot([topStringers.posX],[topStringers.posZ],'.r','markersize',stringerSize) 575 plot([bottomStringers.posX],[bottomStringers.posZ],'.r','markersize',stringerSize) 576 plot([noseTopStringers.posX],[noseTopStringers.posZ],'.r','markersize',stringerSize) 577 plot([noseBottomStringers.posX],[noseBottomStringers.posZ],'.r','markersize',stringerSize) 578 plot(centroid.posX,centroid.posZ,'.k','markerSize',18) 579 plot(sc.posX,sc.posZ,'.g','markersize',18) 580