Exactly what does this script show?

Update on 02 March 2018, 15:34 IST: I have now added another, hopefully better, version of the script (but also kept the old one intact); see in the post below. The new script too comes without comments.


Here is a small little Python script which helps you visualize something about a state of stress in 2D.

If interested in understanding the concept of stress, then do run it, read it, try to understand what it does, and then, if still interested in the concept of stress, try to answer this “simple” little question:

Exactly what does this script show? Exactly what it is that you are visualizing, here?


I had written a few more notes and inline comments in the script, but have deliberately deleted most of them—or at least the ones which might have given you a clue towards answering the above question. I didn’t want to spoil your fun, that’s why.

Once you all finish giving it a try, I will then post another blog-entry here, giving my answer to that question (and in the process, bringing back all the deleted notes and comments).

Anyway, here is the script:


'''
A simple script to help visualize *something* about
a 2D stress tensor.

--Ajit R. Jadhav. Version: 01 March 2018, 21:39 HRS IST.
'''

import math
import numpy as np
import matplotlib.pyplot as plt

# Specifying the input stress
# Note:
# While plotting, we set the x- and y-limits to -150 to +150,
# and enforce the aspect ratio of 1. That is to say, we do not
# allow MatPlotLib to automatically scale the axes, because we
# want to appreciate the changes in the shapes as well sizes in
# the plot.
#
# Therefore, all the input stress-components should be kept
# to within the -100 to +100 (both inclusive) range.
#
# Specify the stress state in this order: xx, xy; yx, yy
# The commas and the semicolon are necessary.

sStress = "-100, 45; 90, 25"

axes = plt.axes()
axes.set_xlim((-150, 150))
axes.set_ylim((-150, 150))
plt.axes().set_aspect('equal', 'datalim')
plt.title(
    "A visualization of *something* about\n" \
    "the 2D stress-state [xx, xy; yx, yy] = [%s]" \
    % sStress)

mStress = np.matrix(sStress)
mStressT = np.transpose(mStress)

mUnitNormal = np.zeros((2, 1))
mTraction = np.zeros((2, 1))

nOrientations = 18
dIncrement = 360.0 / float(nOrientations)
for i in range(0, nOrientations):
    dThetaDegrees = float(i) * dIncrement
    dThetaRads = dThetaDegrees * math.pi / 180.0
    mUnitNormal = [round(math.cos(dThetaRads), 6), round(math.sin(dThetaRads), 6)]
    mTraction = mStressT.dot(mUnitNormal)
    if i == 0:
        plt.plot((0, mTraction[0, 0]), (0, mTraction[0, 1]), 'black', linewidth=1.0)
    else:
        plt.plot((0, mTraction[0, 0]), (0, mTraction[0, 1]), 'gray', linewidth=0.5)
    plt.plot(mTraction[0, 0], mTraction[0, 1], marker='.',
             markeredgecolor='gray', markerfacecolor='gray', markersize=5)
    plt.text(mTraction[0, 0], mTraction[0, 1], '%d' % dThetaDegrees)
    plt.pause(0.05)

plt.show()


Update on 02 March 2018, 15:34 IST:

Here is a second version of a script that does something similar (but continues to lack explanatory comments). One advantage with this version is that you can copy-paste the script to some file, say, MyScript.py, and invoke it from command line, giving the stress components and the number of orientations as command-line inputs, e.g.,

python MyScript.py "100, 0; 0, 50" 12

which makes it easier to try out different states of stress.

The revised code is here:


'''
A simple script to help visualize *something* about
a 2D stress tensor.

--Ajit R. Jadhav. 
History: 
06 March 2018, 10:43 IST: 
In computeTraction(), changed the mUnitNormal code to make it np.matrix() rather than python array
02 March 2018, 15:39 IST; Published the code
'''

import sys
import math
import numpy as np
import matplotlib.pyplot as plt

# Specifying the input stress
# Note:
# While plotting, we set the x- and y-limits to -150 to +150,
# and enforce the aspect ratio of 1. That is to say, we do not
# allow MatPlotLib to automatically scale the axes, because we
# want to appreciate the changes in the shapes as well sizes in
# the plot.
#
# Therefore, all the input stress-components should be kept
# to within the -100 to +100 (both inclusive) range.
#
# Specify the stress state in this order: xx, xy; yx, yy
# The commas and the semicolon are necessary.
# If you run the program from a command-line, you can also
# specify the input stress string in quotes as the first
# command-line argument, and no. of orientations, as the
# second. e.g.:
# python MyScript.py "100, 50; 50, 0" 12
##################################################

gsStress = "-100, 45; 90, 25"
gnOrientations = 18


##################################################

def plotArrow(vTraction, dThetaDegs, clr, axes):
    dx = round(vTraction[0], 6)
    dy = round(vTraction[1], 6)
    if not (math.fabs(dx) < 10e-6 and math.fabs(dy) < 10e-6):
        axes.arrow(0, 0, dx, dy, head_width=3, head_length=9.0, length_includes_head=True, fc=clr, ec=clr)
    axes.annotate(xy=(dx, dy), s='%d' % dThetaDegs, color=clr)


##################################################

def computeTraction(mStressT, dThetaRads):
    vUnitNormal = [round(math.cos(dThetaRads), 6), round(math.sin(dThetaRads), 6)]
    mUnitNormal = np.reshape(vUnitNormal, (2,1))
    mTraction = mStressT.dot(mUnitNormal)
    vTraction = np.squeeze(np.asarray(mTraction))
    return vTraction


##################################################

def main():
    axes = plt.axes()
    axes.set_label("label")
    axes.set_xlim((-150, 150))
    axes.set_ylim((-150, 150))
    axes.set_aspect('equal', 'datalim')
    plt.title(
        "A visualization of *something* about\n" \
        "the 2D stress-state [xx, xy; yx, yy] = [%s]" \
        % gsStress)

    mStress = np.matrix(gsStress)
    mStressT = np.transpose(mStress)
    vTraction = computeTraction(mStressT, 0)
    plotArrow(vTraction, 0, 'red', axes)
    dIncrement = 360.0 / float(gnOrientations)
    for i in range(1, gnOrientations):
        dThetaDegrees = float(i) * dIncrement
        dThetaRads = dThetaDegrees * math.pi / 180.0
        vTraction = computeTraction(mStressT, dThetaRads)
        plotArrow(vTraction, dThetaDegrees, 'gray', axes)
        plt.pause(0.05)
    plt.show()


##################################################

if __name__ == "__main__":
    nArgs = len(sys.argv)
    if nArgs > 1:
        gsStress = sys.argv[1]
    if nArgs > 2:
        gnOrientations = int(sys.argv[2])
    main()



OK, have fun, and if you care to, let me know your answers, guess-works, etc…..


Oh, BTW, I have already taken a version of my last post also to iMechanica, which led to a bit of an interaction there too… However, I had to abruptly cut short all the discussions on the topic because I unexpectedly got way too busy in the affiliation- and accreditation-related work. It was only today that I’ve got a bit of a breather, and so could write this script and this post. Anyway, if you are interested in the concept of stress—issues like what it actually means and all that—then do check out my post at iMechanica, too, here [^].


… Happy Holi, take care to use only safe colors—and also take care not to bother those people who do not want to be bothered by you—by your “play”, esp. the complete strangers…

OK, take care and bye for now. ….


A Song I Like:

(Marathi [Am I right?]) “rang he nave nave…”
Music: Aditya Bedekar
Singer: Shasha Tirupati
Lyrics: Yogesh Damle

 

One thought on “Exactly what does this script show?

  1. Pingback: Stress is defined as the quantity equal to … what? | Ajit Jadhav's Weblog

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

This site uses Akismet to reduce spam. Learn how your comment data is processed.