{"nbformat":4,"nbformat_minor":0,"metadata":{"colab":{"provenance":[],"authorship_tag":"ABX9TyMHSHlAQgOZPGpTSnbIon9F"},"kernelspec":{"name":"python3","display_name":"Python 3"},"language_info":{"name":"python"}},"cells":[{"cell_type":"code","execution_count":31,"metadata":{"id":"BaFgcxSAOffs","executionInfo":{"status":"ok","timestamp":1681511338382,"user_tz":180,"elapsed":5,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}}},"outputs":[],"source":["import numpy as np\n","import time"]},{"cell_type":"code","source":["def triu_1( U, b ):\n","\n","  n = b.shape[ 0 ]\n","  x = np.empty( b.shape )\n","\n","  for i in range( n - 1, -1, -1 ):\n","    tmp = 0.0\n","    for j in range( i + 1, n ):\n","      tmp = tmp + U[ i, j ] * x[ j ]\n","    x[ i ] = ( b[ i ] - tmp ) / U[ i, i ]\n","\n","  return x"],"metadata":{"id":"LA8eEGbDz5Qm","executionInfo":{"status":"ok","timestamp":1681511340687,"user_tz":180,"elapsed":253,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}}},"execution_count":32,"outputs":[]},{"cell_type":"code","source":["n = 1000\n","U = np.random.normal( size = ( n, n ) ).astype( np.float32 )\n","for i in range( n ):\n","  for j in range( i ):\n","    U[ i, j ] = 0.0\n","U[ range( n ), range( n ) ] = U[ range( n ), range( n ) ] + 5.0 * np.sqrt( n )\n","b = np.random.normal( size = ( n, ) ).astype( np.float32 )\n","\n","start = time.time()\n","x = triu_1( U, b )\n","end = time.time()\n","print( end - start )\n","t_1 = end - start\n","\n","print( np.linalg.norm( U @ x - b ) )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"4-9qOWNu1eTT","executionInfo":{"status":"ok","timestamp":1681512262923,"user_tz":180,"elapsed":1848,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"f392a7d5-e4a8-4bf7-8e69-974b2cb5417b"},"execution_count":50,"outputs":[{"output_type":"stream","name":"stdout","text":["1.4638419151306152\n","1.7165360440344537e-14\n"]}]},{"cell_type":"code","source":["print( np.linalg.cond( U ) )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"M27Lqin42MJn","executionInfo":{"status":"ok","timestamp":1681511287618,"user_tz":180,"elapsed":231,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"f8a4d447-1f9f-42a2-c877-1e48387bfd7d"},"execution_count":29,"outputs":[{"output_type":"stream","name":"stdout","text":["81.62477879007007\n"]}]},{"cell_type":"code","source":["def triu_2( U, b ):\n","\n","  n = b.shape[ 0 ]\n","  x = np.empty( b.shape )\n","\n","  for i in range( n - 1, -1, -1 ):\n","    x[ i ] = ( b[ i ] - U[ i, i + 1 : ] @ x[ i + 1 : ] ) / U[ i, i ]\n","\n","  return x"],"metadata":{"id":"7pg1Wssy4VmP","executionInfo":{"status":"ok","timestamp":1681511621544,"user_tz":180,"elapsed":249,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}}},"execution_count":36,"outputs":[]},{"cell_type":"code","source":["start = time.time()\n","x = triu_2( U, b )\n","end = time.time()\n","print( end - start )\n","t_2 = end - start\n","\n","print( np.linalg.norm( U @ x - b ) )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"p-rjFjku5JL3","executionInfo":{"status":"ok","timestamp":1681512271917,"user_tz":180,"elapsed":238,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"206100f9-2064-428f-dea2-514ffcb18fbc"},"execution_count":51,"outputs":[{"output_type":"stream","name":"stdout","text":["0.02025294303894043\n","1.6694613182330644e-14\n"]}]},{"cell_type":"code","source":["print( t_1 / t_2 )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"SPNp24Nn5MXX","executionInfo":{"status":"ok","timestamp":1681512279759,"user_tz":180,"elapsed":2,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"08b73aa3-8b3f-45d6-bb6a-3e3a4e286bae"},"execution_count":52,"outputs":[{"output_type":"stream","name":"stdout","text":["72.27798509658963\n"]}]},{"cell_type":"code","source":["def gauss_solve_1( A, b ):\n","\n","  n = A.shape[ 0 ]\n","\n","  for i in range( n ):\n","    for k in range( i + 1, n ):\n","      alpha = A[ k, i ] / A[ i, i ]\n","      for j in range( i, n ):\n","        A[ k, j ] = A[ k, j ] - alpha * A[ i, j ]\n","      b[ k ] = b[ k ] - alpha * b[ i ]\n","  \n","  x = triu_1( A, b )\n","\n","  return x"],"metadata":{"id":"KdmtPM3y6shH","executionInfo":{"status":"ok","timestamp":1681513607272,"user_tz":180,"elapsed":5,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}}},"execution_count":55,"outputs":[]},{"cell_type":"code","source":["A = np.array(\n","    [\n","        [ 1.0, 1.0, 1.0, 1.0 ],\n","        [ 1.0, 2.0, 3.0, 4.0 ],\n","        [ 1.0, 3.0, 4.0, 3.0 ],\n","        [ 1.0, 1.0, 0.0, 0.0 ]\n","    ]\n",")\n","\n","b = np.array( [ 1.0, 2.0, 3.0, 4.0 ] )\n","\n","x = gauss_solve_1( A.copy(), b.copy() )\n","\n","print( A @ x )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"Wxr5DctvA6ne","executionInfo":{"status":"ok","timestamp":1681513736127,"user_tz":180,"elapsed":230,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"aea49034-6568-496c-cd7a-1bd962edadbb"},"execution_count":56,"outputs":[{"output_type":"stream","name":"stdout","text":["[1. 2. 3. 4.]\n"]}]},{"cell_type":"code","source":["n = 500\n","\n","A = np.random.normal( size = ( n, n ) )\n","b = np.random.normal( size = ( n, ) )\n","\n","start = time.time()\n","x = gauss_solve_1( A.copy(), b.copy() )\n","end = time.time()\n","t_1 = end - start\n","\n","print( t_1 )\n","print( np.linalg.norm( A @ x - b ) )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"kOE-S9ZkB2Sd","executionInfo":{"status":"ok","timestamp":1681515150033,"user_tz":180,"elapsed":25605,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"e94d3be4-7ad9-4910-9578-644a0522d0aa"},"execution_count":84,"outputs":[{"output_type":"stream","name":"stdout","text":["25.328006505966187\n","5.186514203481974e-10\n"]}]},{"cell_type":"code","source":["print( np.linalg.cond( A ) )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"hRdIoHZjCO1m","executionInfo":{"status":"ok","timestamp":1681515162589,"user_tz":180,"elapsed":242,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"544a0b22-84b7-4ee7-e800-b5d4c584ab9d"},"execution_count":85,"outputs":[{"output_type":"stream","name":"stdout","text":["1171.351179029405\n"]}]},{"cell_type":"code","source":["def gauss_solve_2( A, b ):\n","\n","  n = A.shape[ 0 ]\n","\n","  for i in range( n ):\n","    for k in range( i + 1, n ):\n","      alpha = A[ k, i ] / A[ i, i ]\n","      A[ k, i : ] = A[ k, i : ] - alpha * A[ i, i : ]\n","      b[ k ] = b[ k ] - alpha * b[ i ]\n","  \n","  x = triu_2( A, b )\n","\n","  return x"],"metadata":{"id":"BxYpenitCa3U","executionInfo":{"status":"ok","timestamp":1681515165187,"user_tz":180,"elapsed":3,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}}},"execution_count":86,"outputs":[]},{"cell_type":"code","source":["start = time.time()\n","x = gauss_solve_2( A.copy(), b.copy() )\n","end = time.time()\n","t_2 = end - start\n","\n","print( t_2 )\n","print( np.linalg.norm( A @ x - b ) )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"uL2xL6TsD2fM","executionInfo":{"status":"ok","timestamp":1681515167505,"user_tz":180,"elapsed":1105,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"7920b659-85a5-4870-ae67-655f615f84eb"},"execution_count":87,"outputs":[{"output_type":"stream","name":"stdout","text":["0.9684271812438965\n","4.3522229675714344e-10\n"]}]},{"cell_type":"code","source":["def gauss_solve_3( A, b ):\n","\n","  n = A.shape[ 0 ]\n","\n","  for i in range( n ):\n","    alpha = A[ i + 1 :, i ] / A[ i, i ]\n","    A[ i + 1 :, i : ] = A[ i + 1 :, i : ] - alpha.reshape( ( alpha.shape[ 0 ], 1 ) ) * A[ i, i : ]\n","    b[ i + 1 : ] = b[ i + 1 : ] - alpha * b[ i ]\n","  \n","  x = triu_2( A, b )\n","\n","  return x"],"metadata":{"id":"L538_V6gD8OE","executionInfo":{"status":"ok","timestamp":1681515171898,"user_tz":180,"elapsed":228,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}}},"execution_count":88,"outputs":[]},{"cell_type":"code","source":["start = time.time()\n","x = gauss_solve_3( A.copy(), b.copy() )\n","end = time.time()\n","t_3 = end - start\n","\n","print( t_3 )\n","print( np.linalg.norm( A @ x - b ) )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"tdSVZauzEby1","executionInfo":{"status":"ok","timestamp":1681515173921,"user_tz":180,"elapsed":238,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"609862ef-28d6-4fdc-ab67-38fcdf2e66ff"},"execution_count":89,"outputs":[{"output_type":"stream","name":"stdout","text":["0.14435219764709473\n","4.3522229675714344e-10\n"]}]},{"cell_type":"code","source":["print( t_2 / t_3 )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"4eF27hEeFUv8","executionInfo":{"status":"ok","timestamp":1681515175893,"user_tz":180,"elapsed":244,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"0f3d3280-48cd-4870-c683-4259d4cc7fe7"},"execution_count":90,"outputs":[{"output_type":"stream","name":"stdout","text":["6.708780309749495\n"]}]},{"cell_type":"code","source":["print( t_1 / t_3 )"],"metadata":{"colab":{"base_uri":"https://localhost:8080/"},"id":"yUP8v1HQHBg0","executionInfo":{"status":"ok","timestamp":1681515177616,"user_tz":180,"elapsed":373,"user":{"displayName":"Elias Salomão Helou Neto","userId":"05204382143186061711"}},"outputId":"482e2843-9ad4-4062-a51c-8dd923621d7d"},"execution_count":91,"outputs":[{"output_type":"stream","name":"stdout","text":["175.4597915293737\n"]}]},{"cell_type":"code","source":[],"metadata":{"id":"KAnU3VlBH-lb"},"execution_count":null,"outputs":[]}]}