PARADISCYL:How-To-Get-Nodal-Force

Manual 04 for ParaDiS Cylinder Codes

How to get nodal force

In ParaDiS, a dislocation line is represented by piecewise straight segments connecting a set of discrete nodal points. For diagnostic purposes, we may want to find out the driving force on every node. Here we explain how to extract this information from the ParaDiS Cylinder code.

Computing nodal force

Modify makefile in ParaDiS/cylinder directory so that the following line is active.

DEFS += -D_WRITENODEFORCE

Then compile ParaDiS cylinder codes again.

Next, run the executable with the input file, such as that in PARADISCYL:Test-Run,

$ bin/paradiscyl runs/concentric_loop_test.ctrl.

Note that, even if the total simulation step is set as 100 in the control file, the ParaDiS simulation does one time force evaluation and will quit immediately afterwards.

The force data is dumped out as force.out in the outputs/concentric_loop_test/ directory. Each line of this file has 6 values for nodal point coordinate (x, y, z) and nodal force components (f_x, f_y, f_z), respectively.

Because the default stress unit is Pa and length unit is burgMag, the unit of nodal forces is Pa.burgMag^2. Of course, this is the case for cylinder radius R = 1. If we want to interprete the results for the case of R not equal to one, then the force will have to be scaled according to M03 How units are scaled.

---

Plotting nodal force

• 

  Nodal force with all image stress considered

• 

  Nodal force with no BEM image stress

• 

  Nodal force with no Yoffe stress

matlab script

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% plot node force
%
% outdir : directory where nodal force data file is dumped out.
% nodalforcefilename : file name of nodal force data

outdir = '~/Codes/ParaDiS/outputs/concentric_loop_test/';
nodalforcefilename = 'force.out';

fullname_nodalforcefile = sprintf('%s%s',outdir,nodalforcefilename);

data = load(fullname_nodalforcefile);

node_position = data(:,1:3);
nodal_force = data(:,4:6);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% cylinder geometry
r = 1;
theta = linspace(0,2*pi,100);
z = 0;
[x,y,z] = pol2cart(theta,r,z);


figure(1), set(gca,'fontsize',20), hold on
plot(x,y,'-k')                                      % plot cylinder
plot3([node_position(:,1);node_position(1,1)],...   % plot dislocation loop
      [node_position(:,2);node_position(1,2)],...
      [node_position(:,3);node_position(1,3)],'.-')
quiver3(node_position(:,1),node_position(:,2),node_position(:,3),...
        nodal_force(:,1),nodal_force(:,2),nodal_force(:,3),1,'r')
text('position',[0 .07 0],'string','b','fontsize',20,'fontweight','bold')
annotation('arrow', [0.48 0.60],[0.53 0.53],...
           'LineWidth',3,'HeadStyle','plain','HeadWidth',12,'HeadLength',12);
xlabel('x'), ylabel('y')
view(2), axis equal, axis([-1 1 -1 1])
hold off

---

Notes