clear; clear all; clc; format long; E=150; mu=0.4; dispbcelements=[103 0 0 0;113 0 0 0;126 0 0 0;149 0 0 0; 153 0 0 0; 23 0 5 0; 26 0 5 0; 36 0 5 0]; fbcelements(1:22,1:4)=0; for i=1:22 fbcelements(i,1)=i; end fbcelements(24:25,1:4)=0; for i=24:25 fbcelements(i,1)=i; end fbcelements(27:35,1:4)=0; for i=27:35 fbcelements(i,1)=i; end fbcelements(37:102,1:4)=0; for i=37:102 fbcelements(i,1)=i; end fbcelements(104:112,1:4)=0; for i=104:112 fbcelements(i,1)=i; end fbcelements(114:125,1:4)=0; for i=114:125 fbcelements(i,1)=i; end fbcelements(127:148,1:4)=0; for i=127:148 fbcelements(i,1)=i; end fbcelements(150:152,1:4)=0; for i=150:152 fbcelements(i,1)=i; end fbcelements(154:200,1:4)=0; for i=154:200 fbcelements(i,1)=i; end fbcelements(153:153,:)=[]; fbcelements(149:149,:)=[]; fbcelements(126:126,:)=[]; fbcelements(113:113,:)=[]; fbcelements(103:103,:)=[]; fbcelements(36:36,:)=[]; fbcelements(26:26,:)=[]; fbcelements(23:23,:)=[]; xyzofelements=[1 -20.4333 -23.6885 -87.4059; 1 11.3270 -24.0930 -86.0815; 1 -1.3523 -31.8240 -75.2416; 2 -20.8186 -16.3769 -98.9825; 2 -40.1025 -15.3261 -83.5526; 2 -18.0244 -11.4443 -88.4945; 3 -88.0327 -20.2057 124.6373; 3 -49.2881 -16.8016 132.0290; 3 -70.4986 -7.0029 88.3313; 4 -76.7368 -48.0807 61.5649; 4 -38.4633 -48.7990 43.3972; 4 -68.2373 -49.6744 93.8607; 5 -74.5258 -46.5738 143.1554; 5 -93.3111 -45.8081 94.8355; 5 -68.2373 -49.6744 93.8607; 6 -20.8186 -16.3769 -98.9825; 6 -20.4333 -23.6885 -87.4059; 6 -40.1025 -15.3261 -83.5526; 7 -68.2373 -49.6744 93.8607; 7 -43.7603 -49.4117 111.1853; 7 -74.5258 -46.5738 143.1554; 8 -104.5775 -40.1355 122.2658; 8 -74.5258 -46.5738 143.1554; 8 -92.3188 -39.5927 159.8815; 9 -89.3926 -32.1923 160.6916; 9 -106.1549 -31.0810 130.0421; 9 -92.3188 -39.5927 159.8815; 10 -20.4333 -23.6885 -87.4059; 10 -1.3523 -31.8240 -75.2416; 10 -31.1236 -36.0049 -49.5962; 11 -1.3523 -31.8240 -75.2416; 11 0.0810 -38.6876 -52.6897; 11 -31.1236 -36.0049 -49.5962; 12 11.3270 -24.0930 -86.0815; 12 -20.4333 -23.6885 -87.4059; 12 -5.4175 -17.4901 -97.4457; 13 -20.8186 -16.3769 -98.9825; 13 -5.4175 -17.4901 -97.4457; 13 -20.4333 -23.6885 -87.4059; 14 11.3270 -24.0930 -86.0815; 14 26.3618 -33.4085 -59.3037; 14 -1.3523 -31.8240 -75.2416; 15 -5.4175 -17.4901 -97.4457; 15 -18.0244 -11.4443 -88.4945; 15 14.9051 -14.5067 -78.1642; 16 -1.3523 -31.8240 -75.2416; 16 26.3618 -33.4085 -59.3037; 16 0.0810 -38.6876 -52.6897; 17 -56.7616 -24.7525 -57.6934; 17 -20.4333 -23.6885 -87.4059; 17 -31.1236 -36.0049 -49.5962; 18 -77.4892 -3.7099 -32.7977; 18 -98.3911 -11.9581 -11.5482; 18 -97.3285 -3.0388 4.1206; 19 -20.4333 -23.6885 -87.4059; 19 -56.7616 -24.7525 -57.6934; 19 -40.1025 -15.3261 -83.5526; 20 26.3618 -33.4085 -59.3037; 20 11.3270 -24.0930 -86.0815; 20 40.2422 -23.6590 -65.7814; 21 11.3270 -24.0930 -86.0815; 21 -5.4175 -17.4901 -97.4457; 21 14.9051 -14.5067 -78.1642; 22 -20.8186 -16.3769 -98.9825; 22 -18.0244 -11.4443 -88.4945; 22 -5.4175 -17.4901 -97.4457; 23 -45.3464 -44.9662 -3.4297; 23 -34.1323 -48.1268 5.6794; 23 -38.4633 -48.7990 43.3972; 24 -73.8118 -36.3535 -10.4053; 24 -45.3464 -44.9662 -3.4297; 24 -69.2947 -45.5866 25.7985; 25 -16.3498 -44.8691 -14.1788; 25 0.0810 -38.6876 -52.6897; 25 15.4214 -41.3881 -31.4351; 26 -45.3464 -44.9662 -3.4297; 26 -16.3498 -44.8691 -14.1788; 26 -34.1323 -48.1268 5.6794; 27 40.2422 -23.6590 -65.7814; 27 11.3270 -24.0930 -86.0815; 27 14.9051 -14.5067 -78.1642; 28 23.9364 -9.1327 -48.1815; 28 -12.4302 -4.7113 -54.4398; 28 10.6159 -4.2237 -22.0124; 29 -88.0327 -20.2057 124.6373; 29 -98.1258 -7.6650 66.1573; 29 -112.9152 -20.6239 87.2062; 30 -18.0244 -11.4443 -88.4945; 30 -40.1025 -15.3261 -83.5526; 30 -46.7046 -6.5354 -64.5137; 31 0.0810 -38.6876 -52.6897; 31 -16.3498 -44.8691 -14.1788; 31 -31.1236 -36.0049 -49.5962; 32 -56.7616 -24.7525 -57.6934; 32 -31.1236 -36.0049 -49.5962; 32 -73.8118 -36.3535 -10.4053; 33 -31.1236 -36.0049 -49.5962; 33 -16.3498 -44.8691 -14.1788; 33 -45.3464 -44.9662 -3.4297; 34 -40.1025 -15.3261 -83.5526; 34 -56.7616 -24.7525 -57.6934; 34 -71.4192 -13.3841 -51.7763; 35 0.0810 -38.6876 -52.6897; 35 26.3618 -33.4085 -59.3037; 35 15.4214 -41.3881 -31.4351; 36 -16.3498 -44.8691 -14.1788; 36 -38.4633 -48.7990 43.3972; 36 -34.1323 -48.1268 5.6794; 37 -56.7616 -24.7525 -57.6934; 37 -73.8118 -36.3535 -10.4053; 37 -90.4712 -23.3344 -20.7714; 38 -62.0334 -30.1389 170.2612; 38 -24.1336 -26.4878 163.1775; 38 -49.2881 -16.8016 132.0290; 39 -73.8118 -36.3535 -10.4053; 39 -31.1236 -36.0049 -49.5962; 39 -45.3464 -44.9662 -3.4297; 40 -103.3863 -31.6783 15.0978; 40 -73.8118 -36.3535 -10.4053; 40 -96.9559 -41.2047 45.9324; 41 -46.7046 -6.5354 -64.5137; 41 -40.1025 -15.3261 -83.5526; 41 -71.4192 -13.3841 -51.7763; 42 -56.7616 -24.7525 -57.6934; 42 -90.4712 -23.3344 -20.7714; 42 -71.4192 -13.3841 -51.7763; 43 -38.4633 -48.7990 43.3972; 43 -43.7603 -49.4117 111.1853; 43 -68.2373 -49.6744 93.8607; 44 -18.0244 -11.4443 -88.4945; 44 -46.7046 -6.5354 -64.5137; 44 -12.4302 -4.7113 -54.4398; 45 -71.4192 -13.3841 -51.7763; 45 -98.3911 -11.9581 -11.5482; 45 -77.4892 -3.7099 -32.7977; 46 40.2422 -23.6590 -65.7814; 46 57.8010 -30.4582 -31.5163; 46 26.3618 -33.4085 -59.3037; 47 26.3618 -33.4085 -59.3037; 47 57.8010 -30.4582 -31.5163; 47 38.2553 -38.6441 -22.9215; 48 14.9051 -14.5067 -78.1642; 48 23.9364 -9.1327 -48.1815; 48 46.6419 -16.8484 -49.2925; 49 -12.4302 -4.7113 -54.4398; 49 -52.1488 0.5483 -31.8472; 49 -25.6096 0.6197 -2.0068; 50 26.3618 -33.4085 -59.3037; 50 38.2553 -38.6441 -22.9215; 50 15.4214 -41.3881 -31.4351; 51 14.9051 -14.5067 -78.1642; 51 -18.0244 -11.4443 -88.4945; 51 -12.4302 -4.7113 -54.4398; 52 46.6419 -16.8484 -49.2925; 52 40.2422 -23.6590 -65.7814; 52 14.9051 -14.5067 -78.1642; 53 -90.4712 -23.3344 -20.7714; 53 -73.8118 -36.3535 -10.4053; 53 -103.3863 -31.6783 15.0978; 54 23.9364 -9.1327 -48.1815; 54 14.9051 -14.5067 -78.1642; 54 -12.4302 -4.7113 -54.4398; 55 -98.3911 -11.9581 -11.5482; 55 -71.4192 -13.3841 -51.7763; 55 -90.4712 -23.3344 -20.7714; 56 -46.7046 -6.5354 -64.5137; 56 -52.1488 0.5483 -31.8472; 56 -12.4302 -4.7113 -54.4398; 57 -77.4892 -3.7099 -32.7977; 57 -46.7046 -6.5354 -64.5137; 57 -71.4192 -13.3841 -51.7763; 58 40.2422 -23.6590 -65.7814; 58 46.6419 -16.8484 -49.2925; 58 59.1010 -23.2383 -46.0914; 59 -52.1488 0.5483 -31.8472; 59 -46.7046 -6.5354 -64.5137; 59 -77.4892 -3.7099 -32.7977; 60 57.8010 -30.4582 -31.5163; 60 40.2422 -23.6590 -65.7814; 60 59.1010 -23.2383 -46.0914; 61 38.2553 -38.6441 -22.9215; 61 57.8010 -30.4582 -31.5163; 61 55.2373 -37.9617 12.3375; 62 -90.4712 -23.3344 -20.7714; 62 -103.3863 -31.6783 15.0978; 62 -112.7560 -22.3146 19.6431; 63 38.2553 -38.6441 -22.9215; 63 16.1688 -44.4981 18.2126; 63 15.4214 -41.3881 -31.4351; 64 -103.4438 -3.7499 32.4625; 64 -97.3285 -3.0388 4.1206; 64 -111.9131 -12.2067 35.8858; 65 -25.6096 0.6197 -2.0068; 65 -52.1488 0.5483 -31.8472; 65 -61.9642 3.2121 2.7888; 66 -112.7560 -22.3146 19.6431; 66 -103.3863 -31.6783 15.0978; 66 -116.7253 -32.2343 52.6635; 67 -103.3863 -31.6783 15.0978; 67 -96.9559 -41.2047 45.9324; 67 -116.7253 -32.2343 52.6635; 68 -52.1488 0.5483 -31.8472; 68 -77.4892 -3.7099 -32.7977; 68 -78.4442 1.7647 -9.7485; 69 15.4214 -41.3881 -31.4351; 69 16.1688 -44.4981 18.2126; 69 -16.3498 -44.8691 -14.1788; 70 55.2373 -37.9617 12.3375; 70 16.1688 -44.4981 18.2126; 70 38.2553 -38.6441 -22.9215; 71 16.1688 -44.4981 18.2126; 71 55.2373 -37.9617 12.3375; 71 40.3926 -42.2710 58.4336; 72 -98.3911 -11.9581 -11.5482; 72 -90.4712 -23.3344 -20.7714; 72 -112.7560 -22.3146 19.6431; 73 59.1010 -23.2383 -46.0914; 73 69.7321 -22.7243 -13.9361; 73 57.8010 -30.4582 -31.5163; 74 57.8010 -30.4582 -31.5163; 74 69.7321 -22.7243 -13.9361; 74 77.9623 -30.7260 10.3209; 75 69.7321 -22.7243 -13.9361; 75 59.1010 -23.2383 -46.0914; 75 46.6419 -16.8484 -49.2925; 76 -52.1488 0.5483 -31.8472; 76 -78.4442 1.7647 -9.7485; 76 -61.9642 3.2121 2.7888; 77 10.6159 -4.2237 -22.0124; 77 -12.4302 -4.7113 -54.4398; 77 -25.6096 0.6197 -2.0068; 78 -73.8118 -36.3535 -10.4053; 78 -69.2947 -45.5866 25.7985; 78 -96.9559 -41.2047 45.9324; 79 -120.8091 -23.6676 58.1487; 79 -116.7253 -32.2343 52.6635; 79 -118.2629 -31.6563 91.9629; 80 -77.4892 -3.7099 -32.7977; 80 -97.3285 -3.0388 4.1206; 80 -78.4442 1.7647 -9.7485; 81 -45.3464 -44.9662 -3.4297; 81 -38.4633 -48.7990 43.3972; 81 -69.2947 -45.5866 25.7985; 82 -98.3911 -11.9581 -11.5482; 82 -112.7560 -22.3146 19.6431; 82 -111.9131 -12.2067 35.8858; 83 77.9623 -30.7260 10.3209; 83 55.2373 -37.9617 12.3375; 83 57.8010 -30.4582 -31.5163; 84 43.6976 -12.0835 -6.6962; 84 10.6159 -4.2237 -22.0124; 84 11.3158 -4.7213 22.5063; 85 23.9364 -9.1327 -48.1815; 85 43.6976 -12.0835 -6.6962; 85 46.6419 -16.8484 -49.2925; 86 69.7321 -22.7243 -13.9361; 86 46.6419 -16.8484 -49.2925; 86 43.6976 -12.0835 -6.6962; 87 -25.6096 0.6197 -2.0068; 87 -61.9642 3.2121 2.7888; 87 -55.9278 0.5341 43.9214; 88 69.7321 -22.7243 -13.9361; 88 86.8602 -26.1524 18.5976; 88 77.9623 -30.7260 10.3209; 89 -69.2947 -45.5866 25.7985; 89 -76.7368 -48.0807 61.5649; 89 -96.9559 -41.2047 45.9324; 90 23.9364 -9.1327 -48.1815; 90 10.6159 -4.2237 -22.0124; 90 43.6976 -12.0835 -6.6962; 91 -76.7368 -48.0807 61.5649; 91 -68.2373 -49.6744 93.8607; 91 -93.3111 -45.8081 94.8355; 92 39.1237 -13.5858 53.1910; 92 43.6976 -12.0835 -6.6962; 92 11.3158 -4.7213 22.5063; 93 -61.9642 3.2121 2.7888; 93 -78.4442 1.7647 -9.7485; 93 -86.2289 0.7195 27.0344; 94 11.3158 -4.7213 22.5063; 94 10.6159 -4.2237 -22.0124; 94 -25.6096 0.6197 -2.0068; 95 -112.7560 -22.3146 19.6431; 95 -120.8091 -23.6676 58.1487; 95 -111.9131 -12.2067 35.8858; 96 -78.4442 1.7647 -9.7485; 96 -97.3285 -3.0388 4.1206; 96 -86.2289 0.7195 27.0344; 97 -97.3285 -3.0388 4.1206; 97 -98.3911 -11.9581 -11.5482; 97 -111.9131 -12.2067 35.8858; 98 16.1688 -44.4981 18.2126; 98 0.3276 -46.7246 94.9667; 98 -38.4633 -48.7990 43.3972; 99 69.7321 -22.7243 -13.9361; 99 43.6976 -12.0835 -6.6962; 99 72.3139 -21.5623 33.1129; 100 89.2055 -26.4036 35.8411; 100 72.3139 -21.5623 33.1129; 100 79.2496 -28.5255 57.7831; 101 -120.8091 -23.6676 58.1487; 101 -112.7560 -22.3146 19.6431; 101 -116.7253 -32.2343 52.6635; 102 -16.3498 -44.8691 -14.1788; 102 16.1688 -44.4981 18.2126; 102 -38.4633 -48.7990 43.3972; 103 -25.6096 0.6197 -2.0068; 103 -19.3795 -1.5921 45.5447; 103 11.3158 -4.7213 22.5063; 104 -106.1549 -31.0810 130.0421; 104 -89.3926 -32.1923 160.6916; 104 -88.0327 -20.2057 124.6373; 105 86.8602 -26.1524 18.5976; 105 69.7321 -22.7243 -13.9361; 105 72.3139 -21.5623 33.1129; 106 -116.7253 -32.2343 52.6635; 106 -111.6995 -39.3409 83.8489; 106 -118.2629 -31.6563 91.9629; 107 -69.2947 -45.5866 25.7985; 107 -38.4633 -48.7990 43.3972; 107 -76.7368 -48.0807 61.5649; 108 -86.2289 0.7195 27.0344; 108 -98.1258 -7.6650 66.1573; 108 -81.2072 -1.2876 57.3277; 109 -86.2289 0.7195 27.0344; 109 -97.3285 -3.0388 4.1206; 109 -103.4438 -3.7499 32.4625; 110 77.9623 -30.7260 10.3209; 110 86.8602 -26.1524 18.5976; 110 88.0323 -29.6384 34.9591; 111 39.1237 -13.5858 53.1910; 111 4.1889 -8.3091 74.1571; 111 24.1104 -16.9784 100.3874; 112 88.0323 -29.6384 34.9591; 112 86.8602 -26.1524 18.5976; 112 89.2055 -26.4036 35.8411; 113 -19.3795 -1.5921 45.5447; 113 -25.6096 0.6197 -2.0068; 113 -55.9278 0.5341 43.9214; 114 -112.9152 -20.6239 87.2062; 114 -111.9131 -12.2067 35.8858; 114 -120.8091 -23.6676 58.1487; 115 -86.2289 0.7195 27.0344; 115 -55.9278 0.5341 43.9214; 115 -61.9642 3.2121 2.7888; 116 77.9623 -30.7260 10.3209; 116 69.4197 -35.3334 54.0971; 116 55.2373 -37.9617 12.3375; 117 86.8602 -26.1524 18.5976; 117 72.3139 -21.5623 33.1129; 117 89.2055 -26.4036 35.8411; 118 -111.6995 -39.3409 83.8489; 118 -93.3111 -45.8081 94.8355; 118 -104.5775 -40.1355 122.2658; 119 -86.2289 0.7195 27.0344; 119 -81.2072 -1.2876 57.3277; 119 -55.9278 0.5341 43.9214; 120 -103.4438 -3.7499 32.4625; 120 -111.9131 -12.2067 35.8858; 120 -98.1258 -7.6650 66.1573; 121 16.1688 -44.4981 18.2126; 121 40.3926 -42.2710 58.4336; 121 0.3276 -46.7246 94.9667; 122 -70.4986 -7.0029 88.3313; 122 -33.6853 -6.4713 87.4624; 122 -55.9278 0.5341 43.9214; 123 -116.7253 -32.2343 52.6635; 123 -96.9559 -41.2047 45.9324; 123 -111.6995 -39.3409 83.8489; 124 -86.2289 0.7195 27.0344; 124 -103.4438 -3.7499 32.4625; 124 -98.1258 -7.6650 66.1573; 125 -96.9559 -41.2047 45.9324; 125 -93.3111 -45.8081 94.8355; 125 -111.6995 -39.3409 83.8489; 126 11.3158 -4.7213 22.5063; 126 -19.3795 -1.5921 45.5447; 126 4.1889 -8.3091 74.1571; 127 55.2373 -37.9617 12.3375; 127 69.4197 -35.3334 54.0971; 127 40.3926 -42.2710 58.4336; 128 77.9623 -30.7260 10.3209; 128 88.0323 -29.6384 34.9591; 128 69.4197 -35.3334 54.0971; 129 43.6976 -12.0835 -6.6962; 129 39.1237 -13.5858 53.1910; 129 72.3139 -21.5623 33.1129; 130 -38.4633 -48.7990 43.3972; 130 0.3276 -46.7246 94.9667; 130 -43.7603 -49.4117 111.1853; 131 -70.4986 -7.0029 88.3313; 131 -55.9278 0.5341 43.9214; 131 -81.2072 -1.2876 57.3277; 132 72.3139 -21.5623 33.1129; 132 39.1237 -13.5858 53.1910; 132 52.8382 -24.6722 90.2102; 133 40.3926 -42.2710 58.4336; 133 69.4197 -35.3334 54.0971; 133 39.3043 -40.4611 107.6288; 134 -96.9559 -41.2047 45.9324; 134 -76.7368 -48.0807 61.5649; 134 -93.3111 -45.8081 94.8355; 135 -118.2629 -31.6563 91.9629; 135 -111.6995 -39.3409 83.8489; 135 -104.5775 -40.1355 122.2658; 136 -118.2629 -31.6563 91.9629; 136 -112.9152 -20.6239 87.2062; 136 -120.8091 -23.6676 58.1487; 137 89.2055 -26.4036 35.8411; 137 79.2496 -28.5255 57.7831; 137 88.0323 -29.6384 34.9591; 138 -112.9152 -20.6239 87.2062; 138 -98.1258 -7.6650 66.1573; 138 -111.9131 -12.2067 35.8858; 139 -112.9152 -20.6239 87.2062; 139 -106.1549 -31.0810 130.0421; 139 -88.0327 -20.2057 124.6373; 140 69.4197 -35.3334 54.0971; 140 88.0323 -29.6384 34.9591; 140 79.2496 -28.5255 57.7831; 141 40.3926 -42.2710 58.4336; 141 39.3043 -40.4611 107.6288; 141 0.3276 -46.7246 94.9667; 142 0.3276 -46.7246 94.9667; 142 39.3043 -40.4611 107.6288; 142 8.6467 -44.2961 143.4901; 143 4.1889 -8.3091 74.1571; 143 -9.1325 -15.9895 122.0125; 143 24.1104 -16.9784 100.3874; 144 11.3158 -4.7213 22.5063; 144 4.1889 -8.3091 74.1571; 144 39.1237 -13.5858 53.1910; 145 79.2496 -28.5255 57.7831; 145 72.3139 -21.5623 33.1129; 145 52.8382 -24.6722 90.2102; 146 69.4197 -35.3334 54.0971; 146 79.2496 -28.5255 57.7831; 146 61.9612 -33.0296 91.3473; 147 61.9612 -33.0296 91.3473; 147 79.2496 -28.5255 57.7831; 147 52.8382 -24.6722 90.2102; 148 -98.1258 -7.6650 66.1573; 148 -70.4986 -7.0029 88.3313; 148 -81.2072 -1.2876 57.3277; 149 -19.3795 -1.5921 45.5447; 149 -55.9278 0.5341 43.9214; 149 -33.6853 -6.4713 87.4624; 150 -106.1549 -31.0810 130.0421; 150 -112.9152 -20.6239 87.2062; 150 -118.2629 -31.6563 91.9629; 151 61.9612 -33.0296 91.3473; 151 39.3043 -40.4611 107.6288; 151 69.4197 -35.3334 54.0971; 152 -70.4986 -7.0029 88.3313; 152 -98.1258 -7.6650 66.1573; 152 -88.0327 -20.2057 124.6373; 153 -19.3795 -1.5921 45.5447; 153 -33.6853 -6.4713 87.4624; 153 4.1889 -8.3091 74.1571; 154 39.1237 -13.5858 53.1910; 154 24.1104 -16.9784 100.3874; 154 52.8382 -24.6722 90.2102; 155 -70.4986 -7.0029 88.3313; 155 -49.2881 -16.8016 132.0290; 155 -33.6853 -6.4713 87.4624; 156 -43.7603 -49.4117 111.1853; 156 -31.7702 -45.9548 159.0427; 156 -74.5258 -46.5738 143.1554; 157 -35.1892 -35.1717 186.1223; 157 -25.7275 -40.6804 185.6380; 157 -7.2463 -34.9672 177.6640; 158 -74.5258 -46.5738 143.1554; 158 -31.7702 -45.9548 159.0427; 158 -61.3058 -45.6982 175.0607; 159 -106.1549 -31.0810 130.0421; 159 -118.2629 -31.6563 91.9629; 159 -104.5775 -40.1355 122.2658; 160 -43.7603 -49.4117 111.1853; 160 0.3276 -46.7246 94.9667; 160 -31.7702 -45.9548 159.0427; 161 -31.7702 -45.9548 159.0427; 161 0.3276 -46.7246 94.9667; 161 8.6467 -44.2961 143.4901; 162 -31.7702 -45.9548 159.0427; 162 8.6467 -44.2961 143.4901; 162 -0.5949 -41.4291 171.7990; 163 4.1889 -8.3091 74.1571; 163 -33.6853 -6.4713 87.4624; 163 -9.1325 -15.9895 122.0125; 164 61.9612 -33.0296 91.3473; 164 52.8382 -24.6722 90.2102; 164 39.8101 -34.0315 129.9955; 165 8.6467 -44.2961 143.4901; 165 39.3043 -40.4611 107.6288; 165 24.3378 -40.2336 150.3718; 166 61.9612 -33.0296 91.3473; 166 39.8101 -34.0315 129.9955; 166 39.3043 -40.4611 107.6288; 167 -93.3111 -45.8081 94.8355; 167 -74.5258 -46.5738 143.1554; 167 -104.5775 -40.1355 122.2658; 168 -0.5949 -41.4291 171.7990; 168 17.5213 -36.1159 162.5038; 168 -7.2463 -34.9672 177.6640; 169 24.1104 -16.9784 100.3874; 169 -9.1325 -15.9895 122.0125; 169 16.3307 -26.8961 141.9570; 170 -0.5949 -41.4291 171.7990; 170 8.6467 -44.2961 143.4901; 170 24.3378 -40.2336 150.3718; 171 52.8382 -24.6722 90.2102; 171 24.1104 -16.9784 100.3874; 171 16.3307 -26.8961 141.9570; 172 -33.6853 -6.4713 87.4624; 172 -49.2881 -16.8016 132.0290; 172 -9.1325 -15.9895 122.0125; 173 39.3043 -40.4611 107.6288; 173 39.8101 -34.0315 129.9955; 173 24.3378 -40.2336 150.3718; 174 -62.0334 -30.1389 170.2612; 174 -49.2881 -16.8016 132.0290; 174 -88.0327 -20.2057 124.6373; 175 -92.3188 -39.5927 159.8815; 175 -106.1549 -31.0810 130.0421; 175 -104.5775 -40.1355 122.2658; 176 16.3307 -26.8961 141.9570; 176 -9.1325 -15.9895 122.0125; 176 -24.1336 -26.4878 163.1775; 177 -9.1325 -15.9895 122.0125; 177 -49.2881 -16.8016 132.0290; 177 -24.1336 -26.4878 163.1775; 178 -88.0327 -20.2057 124.6373; 178 -89.3926 -32.1923 160.6916; 178 -62.0334 -30.1389 170.2612; 179 -74.5258 -46.5738 143.1554; 179 -75.3256 -43.4635 177.5807; 179 -92.3188 -39.5927 159.8815; 180 39.8101 -34.0315 129.9955; 180 52.8382 -24.6722 90.2102; 180 16.3307 -26.8961 141.9570; 181 39.8101 -34.0315 129.9955; 181 16.3307 -26.8961 141.9570; 181 17.5213 -36.1159 162.5038; 182 -75.3256 -43.4635 177.5807; 182 -54.5417 -40.9243 185.7151; 182 -78.0460 -38.0908 179.8741; 183 -92.3188 -39.5927 159.8815; 183 -75.3256 -43.4635 177.5807; 183 -78.0460 -38.0908 179.8741; 184 -62.0334 -30.1389 170.2612; 184 -89.3926 -32.1923 160.6916; 184 -78.0460 -38.0908 179.8741; 185 -78.0460 -38.0908 179.8741; 185 -89.3926 -32.1923 160.6916; 185 -92.3188 -39.5927 159.8815; 186 -75.3256 -43.4635 177.5807; 186 -61.3058 -45.6982 175.0607; 186 -54.5417 -40.9243 185.7151; 187 17.5213 -36.1159 162.5038; 187 24.3378 -40.2336 150.3718; 187 39.8101 -34.0315 129.9955; 188 -74.5258 -46.5738 143.1554; 188 -61.3058 -45.6982 175.0607; 188 -75.3256 -43.4635 177.5807; 189 17.5213 -36.1159 162.5038; 189 -0.5949 -41.4291 171.7990; 189 24.3378 -40.2336 150.3718; 190 -61.3058 -45.6982 175.0607; 190 -31.7702 -45.9548 159.0427; 190 -54.5417 -40.9243 185.7151; 191 -24.1336 -26.4878 163.1775; 191 -7.2463 -34.9672 177.6640; 191 16.3307 -26.8961 141.9570; 192 -62.0334 -30.1389 170.2612; 192 -35.1892 -35.1717 186.1223; 192 -24.1336 -26.4878 163.1775; 193 -62.0334 -30.1389 170.2612; 193 -78.0460 -38.0908 179.8741; 193 -54.5417 -40.9243 185.7151; 194 16.3307 -26.8961 141.9570; 194 -7.2463 -34.9672 177.6640; 194 17.5213 -36.1159 162.5038; 195 -54.5417 -40.9243 185.7151; 195 -35.1892 -35.1717 186.1223; 195 -62.0334 -30.1389 170.2612; 196 -24.1336 -26.4878 163.1775; 196 -35.1892 -35.1717 186.1223; 196 -7.2463 -34.9672 177.6640; 197 -54.5417 -40.9243 185.7151; 197 -31.7702 -45.9548 159.0427; 197 -25.7275 -40.6804 185.6380; 198 -0.5949 -41.4291 171.7990; 198 -25.7275 -40.6804 185.6380; 198 -31.7702 -45.9548 159.0427; 199 -0.5949 -41.4291 171.7990; 199 -7.2463 -34.9672 177.6640; 199 -25.7275 -40.6804 185.6380; 200 -54.5417 -40.9243 185.7151; 200 -25.7275 -40.6804 185.6380; 200 -35.1892 -35.1717 186.1223]; G=E/(2*(1+mu)); C=1/(16*pi*G*(1-mu)); C1=3-4*mu; C2=1/(8*pi*(1-mu)); C3=1-2*mu; n=2; totalsize=size(xyzofelements); nelement=totalsize(1,1)/3; sizedispbcelements=size(dispbcelements); ndispbcelements=sizedispbcelements(1,1); sizefbcelements=size(fbcelements); nfbcelements=sizefbcelements(1,1); nx(1:nelement)=0; ny(1:nelement)=0; nz(1:nelement)=0; xbar(1:nelement)=0; ybar(1:nelement)=0; zbar(1:nelement)=0; J(1:nelement)=0; for i=1:nelement x1=xyzofelements(3*i-2,2); y1=xyzofelements(3*i-2,3); z1=xyzofelements(3*i-2,4); x2=xyzofelements(3*i-1,2); y2=xyzofelements(3*i-1,3); z2=xyzofelements(3*i-1,4); x3=xyzofelements(3*i,2); y3=xyzofelements(3*i,3); z3=xyzofelements(3*i,4); d=(((y2-y1)*(z3-z1)-(z2-z1)*(y3-y1))^2+((z2-z1)*(x3-x1)-(x2-x1)*(z3-z1))^2+((x2-x1)*(y3-y1)-(y2-y1)*(x3-x1))^2)^0.5; nx(i)=((y2-y1)*(z3-z1)-(z2-z1)*(y3-y1))/d; ny(i)=((z2-z1)*(x3-x1)-(x2-x1)*(z3-z1))/d; nz(i)=((x2-x1)*(y3-y1)-(y2-y1)*(x3-x1))/d; alpha=sqrt((x1-x2)^2+(y1-y2)^2+(z1-z2)^2); beta=sqrt((x2-x3)^2+(y2-y3)^2+(z2-z3)^2); gamma=sqrt((x3-x1)^2+(y3-y1)^2+(z3-z1)^2); sigma=(alpha+beta+gamma)/2; J(i)=2*sqrt(sigma*(sigma-alpha)*(sigma-beta)*(sigma-gamma)); xbar(i)=(x2-x1)*0.25+(x3-x1)*0.5+x1; ybar(i)=(y2-y1)*0.25+(y3-y1)*0.5+y1; zbar(i)=(z2-z1)*0.25+(z3-z1)*0.5+z1; end % For 4 by 4 integration % tonedim=[-0.861136311594 -0.339981043585 0.339981043585 0.861136311594]; % vonedim=[-0.861136311594 -0.339981043585 0.339981043585 0.861136311594]; % weightonedim=[0.347854845137 0.652145154863 0.652145154863 0.347854845137]; % For 8 by 8 integration % tonedim=[-0.960289856498 -0.796666477414 -0.525532409916 -0.183434642496 0.183434642496 0.525532409916 0.796666477414 0.960289856498]; % vonedim=[-0.960289856498 -0.796666477414 -0.525532409916 -0.183434642496 0.183434642496 0.525532409916 0.796666477414 0.960289856498]; % weightonedim=[0.10122853629 0.222381034453 0.313706645878 0.362683783378 0.362683783378 0.313706645878 0.222381034453 0.10122853629]; % 16 by 16 integration tonedim=[-0.989400934992 -0.944575023073 -0.865631202388 -0.755404408355 -0.617876244403 -0.458016777657 -0.281603550779 -0.0950125098376 0.0950125098376 0.281603550779 0.458016777657 0.617876244403 0.755404408355 0.865631202388 0.944575023073 0.989400934992]; vonedim=[-0.989400934992 -0.944575023073 -0.865631202388 -0.755404408355 -0.617876244403 -0.458016777657 -0.281603550779 -0.0950125098376 0.0950125098376 0.281603550779 0.458016777657 0.617876244403 0.755404408355 0.865631202388 0.944575023073 0.989400934992]; weightonedim=[0.027152459411 0.0622535239372 0.0951585116838 0.124628971256 0.149595988817 0.169156519395 0.182603415045 0.189450610455 0.189450610455 0.182603415045 0.169156519395 0.149595988817 0.124628971256 0.0951585116838 0.0622535239372 0.027152459411]; % For 32 by 32 integration % tonedim=[-0.997263861849 -0.985611511545 -0.964762255588 -0.934906075938 -0.896321155766 -0.849367613733 -0.794483795968 -0.73218211874 -0.66304426693 -0.587715757241 -0.506899908932 -0.421351276131 -0.331868602282 -0.239287362252 -0.144471961583 -0.0483076656877 0.0483076656877 0.144471961583 0.239287362252 0.331868602282 0.421351276131 0.506899908932 0.587715757241 0.66304426693 0.73218211874 0.794483795968 0.849367613733 0.896321155766 0.934906075938 0.964762255588 0.985611511545 0.997263861849]; % vonedim=[-0.997263861849 -0.985611511545 -0.964762255588 -0.934906075938 -0.896321155766 -0.849367613733 -0.794483795968 -0.73218211874 -0.66304426693 -0.587715757241 -0.506899908932 -0.421351276131 -0.331868602282 -0.239287362252 -0.144471961583 -0.0483076656877 0.0483076656877 0.144471961583 0.239287362252 0.331868602282 0.421351276131 0.506899908932 0.587715757241 0.66304426693 0.73218211874 0.794483795968 0.849367613733 0.896321155766 0.934906075938 0.964762255588 0.985611511545 0.997263861849]; % weightonedim=[0.00701814576495 0.0162774265831 0.0253910098329 0.0342745478477 0.0428359896785 0.0509978738117 0.0586839394615 0.0658220603578 0.0723456094297 0.078193695762 0.083311711103 0.0876518688047 0.0911736454878 0.0938441590423 0.0956384754512 0.0965398415811 0.0965398415811 0.0956384754512 0.0938441590423 0.0911736454878 0.0876518688047 0.083311711103 0.078193695762 0.0723456094297 0.0658220603578 0.0586839394615 0.0509978738117 0.0428359896785 0.0342745478477 0.0253910098329 0.0162774265831 0.00701814576495]; index=1; t(1:256)=0; v(1:256)=0; weight(1:256)=0; for tindex=1:16 for vindex=1:16 t(index)=tonedim(tindex); v(index)=vonedim(vindex); weight(index)=weightonedim(tindex)*weightonedim(vindex); t(index)=0.5*t(index)+0.5; v(index)=0.5*v(index)+0.5; index=index+1; end end u=t.*(1-v); K=zeros(3*nelement); F(1:3*nelement,1)=0; CPVITxxE(1:nelement)=0; CPVITxyE(1:nelement)=0; CPVITxzE(1:nelement)=0; CPVITyxE(1:nelement)=0; CPVITyyE(1:nelement)=0; CPVITyzE(1:nelement)=0; CPVITzxE(1:nelement)=0; CPVITzyE(1:nelement)=0; CPVITzzE(1:nelement)=0; for i=1:nelement for j=1:nfbcelements x1=xyzofelements(3*fbcelements(j,1)-2,2); y1=xyzofelements(3*fbcelements(j,1)-2,3); z1=xyzofelements(3*fbcelements(j,1)-2,4); x2=xyzofelements(3*fbcelements(j,1)-1,2); y2=xyzofelements(3*fbcelements(j,1)-1,3); z2=xyzofelements(3*fbcelements(j,1)-1,4); x3=xyzofelements(3*fbcelements(j,1),2); y3=xyzofelements(3*fbcelements(j,1),3); z3=xyzofelements(3*fbcelements(j,1),4); x=(x2-x1)*u+(x3-x1)*v+x1; y=(y2-y1)*u+(y3-y1)*v+y1; z=(z2-z1)*u+(z3-z1)*v+z1; ITxxE=0; ITxyE=0; ITxzE=0; ITyxE=0; ITyyE=0; ITyzE=0; ITzxE=0; ITzyE=0; ITzzE=0; IUxxE=0; IUxyE=0; IUxzE=0; IUyxE=0; IUyyE=0; IUyzE=0; IUzxE=0; IUzyE=0; IUzzE=0; for k=1:256 r=sqrt((x(k)-xbar(i))^2+(y(k)-ybar(i))^2+(z(k)-zbar(i))^2); drbydx=(x(k)-xbar(i))/r; drbydy=(y(k)-ybar(i))/r; drbydz=(z(k)-zbar(i))/r; costheta=(1/r)*((x(k)-xbar(i))*nx(fbcelements(j,1))+(y(k)-ybar(i))*ny(fbcelements(j,1))+(z(k)-zbar(i))*nz(fbcelements(j,1))); Txx=(-C2/r^2)*((C3+3*drbydx^2)*costheta); Tyy=(-C2/r^2)*((C3+3*drbydy^2)*costheta); Tzz=(-C2/r^2)*((C3+3*drbydz^2)*costheta); Txy=(-C2/r^2)*(3*drbydx*drbydy*costheta-C3*(ny(fbcelements(j,1))*drbydx-nx(fbcelements(j,1))*drbydy)); Tyx=(-C2/r^2)*(3*drbydy*drbydx*costheta-C3*(nx(fbcelements(j,1))*drbydy-ny(fbcelements(j,1))*drbydx)); Tyz=(-C2/r^2)*(3*drbydy*drbydz*costheta-C3*(nz(fbcelements(j,1))*drbydy-ny(fbcelements(j,1))*drbydz)); Tzy=(-C2/r^2)*(3*drbydz*drbydy*costheta-C3*(ny(fbcelements(j,1))*drbydz-nz(fbcelements(j,1))*drbydy)); Txz=(-C2/r^2)*(3*drbydx*drbydz*costheta-C3*(nz(fbcelements(j,1))*drbydx-nx(fbcelements(j,1))*drbydz)); Tzx=(-C2/r^2)*(3*drbydz*drbydx*costheta-C3*(nx(fbcelements(j,1))*drbydz-nz(fbcelements(j,1))*drbydx)); ITxxE=ITxxE+Txx*(1-v(k))*J(fbcelements(j,1))*weight(k); ITxyE=ITxyE+Txy*(1-v(k))*J(fbcelements(j,1))*weight(k); ITxzE=ITxzE+Txz*(1-v(k))*J(fbcelements(j,1))*weight(k); ITyxE=ITyxE+Tyx*(1-v(k))*J(fbcelements(j,1))*weight(k); ITyyE=ITyyE+Tyy*(1-v(k))*J(fbcelements(j,1))*weight(k); ITyzE=ITyzE+Tyz*(1-v(k))*J(fbcelements(j,1))*weight(k); ITzxE=ITzxE+Tzx*(1-v(k))*J(fbcelements(j,1))*weight(k); ITzyE=ITzyE+Tzy*(1-v(k))*J(fbcelements(j,1))*weight(k); ITzzE=ITzzE+Tzz*(1-v(k))*J(fbcelements(j,1))*weight(k); Uxx=(C/r)*(C1+drbydx^2); Uyy=(C/r)*(C1+drbydy^2); Uzz=(C/r)*(C1+drbydz^2); Uxy=(C/r)*drbydx*drbydy; Uyx=Uxy; Uyz=(C/r)*drbydy*drbydz; Uzy=Uyz; Uxz=(C/r)*drbydx*drbydz; Uzx=Uxz; IUxxE=IUxxE+Uxx*(1-v(k))*J(fbcelements(j,1))*weight(k); IUxyE=IUxyE+Uxy*(1-v(k))*J(fbcelements(j,1))*weight(k); IUxzE=IUxzE+Uxz*(1-v(k))*J(fbcelements(j,1))*weight(k); IUyxE=IUyxE+Uyx*(1-v(k))*J(fbcelements(j,1))*weight(k); IUyyE=IUyyE+Uyy*(1-v(k))*J(fbcelements(j,1))*weight(k); IUyzE=IUyzE+Uyz*(1-v(k))*J(fbcelements(j,1))*weight(k); IUzxE=IUzxE+Uzx*(1-v(k))*J(fbcelements(j,1))*weight(k); IUzyE=IUzyE+Uzy*(1-v(k))*J(fbcelements(j,1))*weight(k); IUzzE=IUzzE+Uzz*(1-v(k))*J(fbcelements(j,1))*weight(k); end ITxxE=0.5*ITxxE; ITxyE=0.5*ITxyE; ITxzE=0.5*ITxzE; ITyxE=0.5*ITyxE; ITyyE=0.5*ITyyE; ITyzE=0.5*ITyzE; ITzxE=0.5*ITzxE; ITzyE=0.5*ITzyE; ITzzE=0.5*ITzzE; IUxxE=0.5*IUxxE; IUxyE=0.5*IUxyE; IUxzE=0.5*IUxzE; IUyxE=0.5*IUyxE; IUyyE=0.5*IUyyE; IUyzE=0.5*IUyzE; IUzxE=0.5*IUzxE; IUzyE=0.5*IUzyE; IUzzE=0.5*IUzzE; CPVITxxE(i)=CPVITxxE(i)-ITxxE; CPVITxyE(i)=CPVITxyE(i)-ITxyE; CPVITxzE(i)=CPVITxzE(i)-ITxzE; CPVITyxE(i)=CPVITyxE(i)-ITyxE; CPVITyyE(i)=CPVITyyE(i)-ITyyE; CPVITyzE(i)=CPVITyzE(i)-ITyzE; CPVITzxE(i)=CPVITzxE(i)-ITzxE; CPVITzyE(i)=CPVITzyE(i)-ITzyE; CPVITzzE(i)=CPVITzzE(i)-ITzzE; if i==fbcelements(j,1) CPVITxxE(i)=CPVITxxE(i)+ITxxE; CPVITxyE(i)=CPVITxyE(i)+ITxyE; CPVITxzE(i)=CPVITxzE(i)+ITxzE; CPVITyxE(i)=CPVITyxE(i)+ITyxE; CPVITyyE(i)=CPVITyyE(i)+ITyyE; CPVITyzE(i)=CPVITyzE(i)+ITyzE; CPVITzxE(i)=CPVITzxE(i)+ITzxE; CPVITzyE(i)=CPVITzyE(i)+ITzyE; CPVITzzE(i)=CPVITzzE(i)+ITzzE; end K(3*i-2,3*fbcelements(j,1)-2)=ITxxE; K(3*i-2,3*fbcelements(j,1)-1)=ITxyE; K(3*i-2,3*fbcelements(j,1))=ITxzE; K(3*i-1,3*fbcelements(j,1)-2)=ITyxE; K(3*i-1,3*fbcelements(j,1)-1)=ITyyE; K(3*i-1,3*fbcelements(j,1))=ITyzE; K(3*i,3*fbcelements(j,1)-2)=ITzxE; K(3*i,3*fbcelements(j,1)-1)=ITzyE; K(3*i,3*fbcelements(j,1))=ITzzE; F(3*i-2)=F(3*i-2)+fbcelements(j,2)*IUxxE+fbcelements(j,3)*IUxyE+fbcelements(j,4)*IUxzE; F(3*i-1)=F(3*i-1)+fbcelements(j,2)*IUyxE+fbcelements(j,3)*IUyyE+fbcelements(j,4)*IUyzE; F(3*i)=F(3*i)+fbcelements(j,2)*IUzxE+fbcelements(j,3)*IUzyE+fbcelements(j,4)*IUzzE; end for j=1:ndispbcelements x1=xyzofelements(3*dispbcelements(j,1)-2,2); y1=xyzofelements(3*dispbcelements(j,1)-2,3); z1=xyzofelements(3*dispbcelements(j,1)-2,4); x2=xyzofelements(3*dispbcelements(j,1)-1,2); y2=xyzofelements(3*dispbcelements(j,1)-1,3); z2=xyzofelements(3*dispbcelements(j,1)-1,4); x3=xyzofelements(3*dispbcelements(j,1),2); y3=xyzofelements(3*dispbcelements(j,1),3); z3=xyzofelements(3*dispbcelements(j,1),4); x=(x2-x1)*u+(x3-x1)*v+x1; y=(y2-y1)*u+(y3-y1)*v+y1; z=(z2-z1)*u+(z3-z1)*v+z1; ITxxE=0; ITxyE=0; ITxzE=0; ITyxE=0; ITyyE=0; ITyzE=0; ITzxE=0; ITzyE=0; ITzzE=0; IUxxE=0; IUxyE=0; IUxzE=0; IUyxE=0; IUyyE=0; IUyzE=0; IUzxE=0; IUzyE=0; IUzzE=0; for k=1:256 r=sqrt((x(k)-xbar(i))^2+(y(k)-ybar(i))^2+(z(k)-zbar(i))^2); drbydx=(x(k)-xbar(i))/r; drbydy=(y(k)-ybar(i))/r; drbydz=(z(k)-zbar(i))/r; costheta=(1/r)*((x(k)-xbar(i))*nx(dispbcelements(j,1))+(y(k)-ybar(i))*ny(dispbcelements(j,1))+(z(k)-zbar(i))*nz(dispbcelements(j,1))); Txx=(-C2/r^2)*((C3+3*drbydx^2)*costheta); Tyy=(-C2/r^2)*((C3+3*drbydy^2)*costheta); Tzz=(-C2/r^2)*((C3+3*drbydz^2)*costheta); Txy=(-C2/r^2)*(3*drbydx*drbydy*costheta-C3*(ny(dispbcelements(j,1))*drbydx-nx(dispbcelements(j,1))*drbydy)); Tyx=(-C2/r^2)*(3*drbydy*drbydx*costheta-C3*(nx(dispbcelements(j,1))*drbydy-ny(dispbcelements(j,1))*drbydx)); Tyz=(-C2/r^2)*(3*drbydy*drbydz*costheta-C3*(nz(dispbcelements(j,1))*drbydy-ny(dispbcelements(j,1))*drbydz)); Tzy=(-C2/r^2)*(3*drbydz*drbydy*costheta-C3*(ny(dispbcelements(j,1))*drbydz-nz(dispbcelements(j,1))*drbydy)); Txz=(-C2/r^2)*(3*drbydx*drbydz*costheta-C3*(nz(dispbcelements(j,1))*drbydx-nx(dispbcelements(j,1))*drbydz)); Tzx=(-C2/r^2)*(3*drbydz*drbydx*costheta-C3*(nx(dispbcelements(j,1))*drbydz-nz(dispbcelements(j,1))*drbydx)); ITxxE=ITxxE+Txx*(1-v(k))*J(dispbcelements(j,1))*weight(k); ITxyE=ITxyE+Txy*(1-v(k))*J(dispbcelements(j,1))*weight(k); ITxzE=ITxzE+Txz*(1-v(k))*J(dispbcelements(j,1))*weight(k); ITyxE=ITyxE+Tyx*(1-v(k))*J(dispbcelements(j,1))*weight(k); ITyyE=ITyyE+Tyy*(1-v(k))*J(dispbcelements(j,1))*weight(k); ITyzE=ITyzE+Tyz*(1-v(k))*J(dispbcelements(j,1))*weight(k); ITzxE=ITzxE+Tzx*(1-v(k))*J(dispbcelements(j,1))*weight(k); ITzyE=ITzyE+Tzy*(1-v(k))*J(dispbcelements(j,1))*weight(k); ITzzE=ITzzE+Tzz*(1-v(k))*J(dispbcelements(j,1))*weight(k); Uxx=(C/r)*(C1+drbydx^2); Uyy=(C/r)*(C1+drbydy^2); Uzz=(C/r)*(C1+drbydz^2); Uxy=(C/r)*drbydx*drbydy; Uyx=Uxy; Uyz=(C/r)*drbydy*drbydz; Uzy=Uyz; Uxz=(C/r)*drbydx*drbydz; Uzx=Uxz; IUxxE=IUxxE+Uxx*(1-v(k))*J(dispbcelements(j,1))*weight(k); IUxyE=IUxyE+Uxy*(1-v(k))*J(dispbcelements(j,1))*weight(k); IUxzE=IUxzE+Uxz*(1-v(k))*J(dispbcelements(j,1))*weight(k); IUyxE=IUyxE+Uyx*(1-v(k))*J(dispbcelements(j,1))*weight(k); IUyyE=IUyyE+Uyy*(1-v(k))*J(dispbcelements(j,1))*weight(k); IUyzE=IUyzE+Uyz*(1-v(k))*J(dispbcelements(j,1))*weight(k); IUzxE=IUzxE+Uzx*(1-v(k))*J(dispbcelements(j,1))*weight(k); IUzyE=IUzyE+Uzy*(1-v(k))*J(dispbcelements(j,1))*weight(k); IUzzE=IUzzE+Uzz*(1-v(k))*J(dispbcelements(j,1))*weight(k); end ITxxE=0.5*ITxxE; ITxyE=0.5*ITxyE; ITxzE=0.5*ITxzE; ITyxE=0.5*ITyxE; ITyyE=0.5*ITyyE; ITyzE=0.5*ITyzE; ITzxE=0.5*ITzxE; ITzyE=0.5*ITzyE; ITzzE=0.5*ITzzE; IUxxE=0.5*IUxxE; IUxyE=0.5*IUxyE; IUxzE=0.5*IUxzE; IUyxE=0.5*IUyxE; IUyyE=0.5*IUyyE; IUyzE=0.5*IUyzE; IUzxE=0.5*IUzxE; IUzyE=0.5*IUzyE; IUzzE=0.5*IUzzE; CPVITxxE(i)=CPVITxxE(i)-ITxxE; CPVITxyE(i)=CPVITxyE(i)-ITxyE; CPVITxzE(i)=CPVITxzE(i)-ITxzE; CPVITyxE(i)=CPVITyxE(i)-ITyxE; CPVITyyE(i)=CPVITyyE(i)-ITyyE; CPVITyzE(i)=CPVITyzE(i)-ITyzE; CPVITzxE(i)=CPVITzxE(i)-ITzxE; CPVITzyE(i)=CPVITzyE(i)-ITzyE; CPVITzzE(i)=CPVITzzE(i)-ITzzE; if i==fbcelements(j,1) CPVITxxE(i)=CPVITxxE(i)+ITxxE; CPVITxyE(i)=CPVITxyE(i)+ITxyE; CPVITxzE(i)=CPVITxzE(i)+ITxzE; CPVITyxE(i)=CPVITyxE(i)+ITyxE; CPVITyyE(i)=CPVITyyE(i)+ITyyE; CPVITyzE(i)=CPVITyzE(i)+ITyzE; CPVITzxE(i)=CPVITzxE(i)+ITzxE; CPVITzyE(i)=CPVITzyE(i)+ITzyE; CPVITzzE(i)=CPVITzzE(i)+ITzzE; end K(3*i-2,3*dispbcelements(j,1)-2)=IUxxE; K(3*i-2,3*dispbcelements(j,1)-1)=IUxyE; K(3*i-2,3*dispbcelements(j,1))=IUxzE; K(3*i-1,3*dispbcelements(j,1)-2)=IUyxE; K(3*i-1,3*dispbcelements(j,1)-1)=IUyyE; K(3*i-1,3*dispbcelements(j,1))=IUyzE; K(3*i,3*dispbcelements(j,1)-2)=IUzxE; K(3*i,3*dispbcelements(j,1)-1)=IUzyE; K(3*i,3*dispbcelements(j,1))=IUzzE; F(3*i-2)=F(3*i-2)+dispbcelements(j,2)*ITxxE+dispbcelements(j,3)*ITxyE+dispbcelements(j,4)*ITxzE; F(3*i-1)=F(3*i-1)+dispbcelements(j,2)*ITyxE+dispbcelements(j,3)*ITyyE+dispbcelements(j,4)*ITyzE; F(3*i)=F(3*i)+dispbcelements(j,2)*ITzxE+dispbcelements(j,3)*ITzyE+dispbcelements(j,4)*ITzzE; if i==dispbcelements(j,1) F(3*i-2)=F(3*i-2)-dispbcelements(j,2)*ITxxE-dispbcelements(j,3)*ITxyE-dispbcelements(j,4)*ITxzE; F(3*i-1)=F(3*i-1)-dispbcelements(j,2)*ITyxE-dispbcelements(j,3)*ITyyE-dispbcelements(j,4)*ITyzE; F(3*i)=F(3*i)-dispbcelements(j,2)*ITzxE-dispbcelements(j,3)*ITzyE-dispbcelements(j,4)*ITzzE; end end end for i=1:nelement for j=1:nfbcelements if i==fbcelements(j,1) K(3*i-2,3*fbcelements(j,1)-2)=CPVITxxE(i); K(3*i-2,3*fbcelements(j,1)-1)=CPVITxyE(i); K(3*i-2,3*fbcelements(j,1))=CPVITxzE(i); K(3*i-1,3*fbcelements(j,1)-2)=CPVITyxE(i); K(3*i-1,3*fbcelements(j,1)-1)=CPVITyyE(i); K(3*i-1,3*fbcelements(j,1))=CPVITyzE(i); K(3*i,3*fbcelements(j,1)-2)=CPVITzxE(i); K(3*i,3*fbcelements(j,1)-1)=CPVITzyE(i); K(3*i,3*fbcelements(j,1))=CPVITzzE(i); end end for j=1:ndispbcelements if i==dispbcelements(j,1) F(3*i-2)=F(3*i-2)+dispbcelements(j,2)*CPVITxxE(i)+dispbcelements(j,3)*CPVITxyE(i)+dispbcelements(j,4)*CPVITxzE(i); F(3*i-1)=F(3*i-1)+dispbcelements(j,2)*CPVITyxE(i)+dispbcelements(j,3)*CPVITyyE(i)+dispbcelements(j,4)*CPVITyzE(i); F(3*i)=F(3*i)+dispbcelements(j,2)*CPVITzxE(i)+dispbcelements(j,3)*CPVITzyE(i)+dispbcelements(j,4)*CPVITzzE(i); end end end %Solve U=K\F; % Display unknown dofs for each element U