Why is my neural net getting a poor accuracy? EPOCH=40
batch_size=50
mnist=input_data.read_data_sets("MNIST_data", one_hot=True)

tf.reset_default_graph()

input_X = tf.placeholder(tf.float32, shape=[None, 784])
input_y = tf.placeholder(tf.int64, shape=[None, 10])

input_layer=tf.layers.dense(input_X, 784, activation=tf.nn.sigmoid)
hidden1=tf.layers.dense(input_layer, 256, activation=tf.nn.sigmoid)
hidden2=tf.layers.dense(hidden1, 256, activation=tf.nn.sigmoid)
output=tf.layers.dense(hidden2, units=10)
output=tf.nn.softmax(output)


entropy=tf.nn.softmax_cross_entropy_with_logits(labels=input_y, logits=output)
loss=tf.reduce_mean(entropy)
step=tf.train.GradientDescentOptimizer(0.01).minimize(loss)

#correct=tf.nn.in_top_k(tf.argmax(output, y, 1) rank error

correct=tf.equal(tf.argmax(output, 1), tf.argmax(input_y, 1))
accuracy = tf.reduce_mean(tf.cast(correct, tf.float32))


sess=tf.InteractiveSession()
sess.run(tf.global_variables_initializer())

for epoch in range(EPOCH):
  for i in range(1000):
    X, y=mnist.train.next_batch(batch_size)
    sess.run(loss, feed_dict={input_X:X, input_y:y})
  acc_train=accuracy.eval(feed_dict={input_X:X, input_y:y})
  acc_val=accuracy.eval(feed_dict={input_X: mnist.validation.images
                                  , input_y: mnist.validation.labels})
  print(epoch, "Train accuracy: ", acc_train, "\n Val accuracy: ", acc_val)


Results went something like this-

0 Train accuracy:  0.1 
 Val accuracy:  0.0986
1 Train accuracy:  0.14 
 Val accuracy:  0.0986
2 Train accuracy:  0.06 
 Val accuracy:  0.0986
3 Train accuracy:  0.06 
 Val accuracy:  0.0986


I have referenced part of the code from O'Reily textbook on Tensorflow.

EDIT- found this O' Reily guide online-https://github.com/ageron/handson-ml/blob/master/10_introduction_to_artificial_neural_networks.ipynb

They have a code similar to mine but their's works just fine.
Reason is you did not train the network.

You need to run the step op in your code. You defined it, but didn't use it anywhere later. Running the  loss only calculates the loss, but does not train your network.

after fixing your code (remove duplicate softmax, using adadelta, run step) you should get: